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Chapter  1 
INTRODUCTION 

BACKGROUND 

This  report  describes  the  modification  of  the  Vulnerability  Model 
(VM)  ,  a  computer  simulation  which  provides  quantitative  measures  of  the 
consequences  of  maritime  spills  of  hazardous  materials.  The  VM  is  being 
developed  for  the  U.S.  Coast  Guard  under  contract  DOT-CG-33377-A.  The 
first,  second,  and  third  stages  of  development  are  described  in  refer¬ 
ences  (1J,  (2),  and  (31,  respectively.  This  present  report  describes 
the  latest  modifications. 

The  VM  is  a  research  tool,  one  use  for  which  is  in  the  USCG  Risk 
Management  Program.  It  has  been  designed  to  treat  virtually  all  of  the 
large  class  of  hazardous  materials  carried  in  bulk  in  marine  transport. 
The  simulation  starts  with  a  description  of  the  nature  of  the  spill, 
simulates  the  dispersion  of  the  hazardous  material,  and  assesses  the 
immediate  effects  of  the  spill  on  surrounding  vulnerable  resources, 
namely:  people,  property,  and  the  environment. 


The  VM  requires  three  types  of  descriptive  data  that  define:  (1)  the 
spill,  (2)  the  physical  setting  in  whicn  the  spill  occurs,  and  (3)  the 
vulnerable  resources  that  are  subject  to  the  effects  of  the  spill.  The 
spill  is  described  in  terms  of  its  location  and  spill  rate,  the  physical 
and  chemical  properties  of  the  spilled  material,  and  the  quantity  of  the 
spill.  The  physical  setting  is  described  in  terms  of  the  geometric  con¬ 
figuration  of  the  shoreline (s) ,  hydrologic/oceanographic  properties,  and 
meteorological  data.  Vulnerable  resources  are  described  in  terms  of 


(1)  Eisenberg,  N.  A.,  C.  J.  Lynch,  and  R.  J.  Breeding,  Vulnerability 
Model:  A  Simulation  System  for  Assessing  Damage  Resulting  from 
Marine  Spills,  CG-D-136-75  (NTIS  AD-A015  245),  Final  Report,  pre¬ 
pared  by  Enviro  Control,  Inc.,  for  Department  of  Transportation, 

U.S.  Coast  Guard,  Office  of  Research  and  Development,  June  1975. 

(2)  Rausch,  A.  H. ,  N.  A.  Eisenberg,  and  C.  J.  Lynch,  Continuing  Develop¬ 
ment  of  the  Vulnerability  Model:  A  Simulation  System  for  Assessing 
Damage  Resulting  from  Marine  Spills,  Final  Report,  prepared  by 
Enviro  Control,  Inc.,  for  Department  of  Transportation,  U.S.  Coast 
Guard,  Office  of  Research  and  Development,  February  1977. 

(3)  Rausch,  A.  H.,  C.  K.  Tsao,  and  R.  M.  Rowley,  Third-Stage  Develop¬ 
ment  of  the  Vulnerability  Model:  A  Simulation  System  for  Assessing 
Damage  Resulting  from  Marine  Spills,  Final  Report,  prepared  by 
Enviro  Control,  Inc.,  for  Department  of  Transportation,  U.S.  Coast 
Guard,  Office  of  Research  and  Development,  June  1977. 
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demographic  distribution,  property  distribution,  and  land/water  use. 
The  geographic  area  of  concern  may  represent  any  user-defined  location. 
The  physical  setting  and  the  distribution  of  vulnerable  resources  are 
described  in  terms  of  mutually  exclusive  geographic  cells  that  cover 
the  entire  area  of  concern. 


The  VM  operates  in  two  phases.  Phase  I  simulates  thu  spill,  the 
physical  and  chemical  transformations  of  the  spilled  substance,  and  its 
dissemination  in  space.  This  phase  covers  the  time  period  from  the 
initiation  of  the  spill  until  a  user-specified  time  has  elapsed .  Phase 
I  consists  of  submodels  interconnected  by  an  executive  routine,  with 
built-in  logic  dictating  the  sequence  of  submodel  processing  as  a  func¬ 
tion  of  tho  spill  development.  Submodels  depicting  spill  development 
simulate  the  following  phenomena!  (1)  cargo  v»  .ng,  (2)  surface 
spreading  (with  or  without  evaporation),  (3)  water  mixing,  (4)  sinking 
and  boiling,  (5)  air  dispersion,  and  (6)  fire  and  explosion.  A  time- 
history  file  of  the  spill  sequence  simulated  during  the  first  phase  is 
retained  in  computer  storage  on  magnetic  tape  and  disk. 


In  Phase  II  the  computer  first  matches  this  time-history  file  to 
the  vulnerable  resources  map,  and  then  assesses  the  effects  of  toxicity, 
explosion  and/or  fire  on  the  vulnerable  resources  as  a  function  of  time. 
Estimates  of  deaths  and  nonlethal  injuries  to  people  and  of  damage  to 
property  are  presented  in  computer-generated  tables.  A  summary  of  the 
types  of  Phase  II  damage  is  given  in  the  following  table. 


PHASE  II  DAMAGE  ASSESSMENT 


DAMAGE-CAUSING 

EVENT 

VULNERABLE 

RESOURCE 

TYPE  or 

INJURY  OR  DAMAGE 

cause  or 

INJURY  OR  DAMAGE 

TOXICITY 

People 

Death 

Nonlethal  Injury 

Irritation 

Toxic  vapor i 

Concentration  or 
cumulative  doet> 

EXPLOSION 

people 

Death 

Honiethai  injury 

•  Eardrum  rupture 

•  Dona  fracture 

•  Puncture  wound 

•  Multiple  injury 

Direct  Blaet i  impact 

Direct  Blaet 

Impact 

Plying  Fragment! 

Two  or  more  of  the  above 

Structural 

Structural  Damage 

Giaia  Breakage 

Direct  Blaet 

POOL  BURNING 
i  FIREBALL 

FLASH  FIRE 

People 

Death 

Flret-Degree  Burn 

Thermal  Radiation 

Structural 

ignition 
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SCOPE  OF  WORK 


The  work  on  the  Vulnerability  Model  described  in  this  report  is  made 
up  of  three  tasks. 

Task  1 — Modification  of  plume  and  puff  models  for  toxic  damage 
assessment. 

A.  Remodeling  of  the  plume  model  to  reduce  computing 
time  and  storage  space. 

0.  Computing  the  dosage  for  the  puff  model  by  direct 
integration. 

C.  Transferring  dosage  values  computed  in  Phase  I  as 
specified  in  (A)  and  (B)  above  directly  to  Phase  II. 

Task  2 — Modification  of  the  thermal  injury  and  lethality  criteria 
in  Phase  II. 

A.  Rederivation  of  the  dose  equation  based  on  applicable 
experimental  data. 

B.  Revision  of  the  dose  criteria  for  thermal  injury 
and  lethality. 

C.  Recomputation  of  the  coefficients  in  the  probit  equation. 

Task  3 — Remodeling  of  flash  fire. 

A.  Remodeling  of  the  puff  model  using  an  adaptation  of 
the  fireball  model. 

B.  Remodeling  of  the  plume  model  on  the  same  basis  as 
in  (A)  above. 

C.  Recalculation  of  the  mass  burned  in  flash  fire. 

The  results  of  the  three  tasks  are  presented  in  order  in  the  follow¬ 
ing  three  chapters  (Chapters  2,  3,  and  4).  The  final  chapter  (5)  sum¬ 
marizes  the  revisions  to  the  VM  computer  program  incorporating  the  modi¬ 
fications  that  have  been  developed. 


CONCLUSIONS 

The  modifications  described  in  this  report  make  concrete  improve¬ 
ments  to  the  VM.  The  direct  integration  of  toxic  dose  improves  the  ac¬ 
curacy  of  results  and  saves  computing  time  and  data  storage.  The  modi¬ 
fication  of  the  probit  equation  for  thermal  injury  and  lethality  makes 
the  calculation  consistent  with  experimental  results.  The  modification 
of  the  flash  fire  model  provides  a  more  accurate  simulation  of  the  ther¬ 
mal  radiation  emitted  from  burning  vapor  clouds. 
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Chapter  2 

MODIFICATION  OF  PLUME  AND  PUFF  MODELS 


INTRODUCTION 

the  duration  of 
The  general  form  of 


(2-1) 


number.  To  evaluate 
been  used  in  the 
vapor  concentration 
average  of  the  nth 
power  of  concentrations  over  time.  The  drawbacks  of  the  finite  differ¬ 
ence  model  are: 

•  Very  few  physical  inferences  can  be  drawn  from  the  numerical 
results. 

•  The  calculation  must  start  from  the  first  time  step  and 
follow  the  time-step  sequence.  There  is  no  way  to  calculate 
the  concentration  at  the  (j  +1)  time  step  without  calculating 
the  concentration  of  the  first  j  time  steps. 

•  The  accuracy  of  the  results  is  primarily  dependent  upon  the 
si2e  of  the  mesh  and  the  fineness  of  the  time  step. 

•  It  prolongs  the  computing  time  and  burdens  the  storage 
capacity  of  the  computer. 

To  improve  accuracy  and  to  reduce  computing  time,  Equation  (2-1)  is 
now  calculated  by  direct  integration  rather  than  time-stepped  integra¬ 
tion.  This  chapter  describes  how  the  direct  integration  is  accomplished 
in  the  plume  and  puff  models  for  both  exterior  and  interior  dosages. 


For  toxic  injury,  the  dose  v  depends  upon  both 
exposure  and  the  concentration  level  experienced, 
the  dose  is  as  follows: 


v  = 


Cn  dt 


where  C  is  the  vapor  concentration  and  n  is  a  real 
the  dose,  a  finite  difference  model  had  previously 
Vulnerability  Model  (VM) .  It  first  calculated  the 
in  each  cell  at  each  time  step  and  then  summed  the 


V 
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DERIVATION  OF  EXTERIOR  DOSAGE 

In  the  VM  the  simulation  of  vapor  dispersion  is  based  on  the  diffu¬ 
sion  equation  (4,5): 


3£  3C  _  J_ 

at  u  ax  "  ax 


D  5£l+  — 

3 

i  —  - 

L  li 

x  3x]  3y 

D 

[  y  3yJ 

3z 

Vz  3z. 

(2-2) 


where  U  is  tne  air  current  or  wind  velocity  along  the  x-axis,  and  Dx  and 
Dy  and  Dz  ace  the  diffusion  coefficients.  Consider  an  instantaneous 
point  source  of  strength  m  located  at  the  origin,  with  the  boundary  on 
the  ground  (z  =0),  the  wind  velocity  V  and  the  diffusion  coefficients 
constant;  then  the  solution  of  Equation  (2-2)  is: 


C 


i 


_ 2m _ 

(2ir)3/2oxay  o2 


exp 


(x  -Ut)2 

2  ax2 


(2-3) 


where  ox  =  /2  Dxt  and  Gy  =  /2  Dyt  and  oz  =  \/2  Dzt  are  the  dispersion 
coefficients,  and 


m  =  Iff  — dxdydz 

-CO  2 


(2- 3a) 


is  the  total  amount  of  material  released  at  the  origin  at  time  (t  =0). 
Equation  (2-3)  is  the  puff  solution. 

For  a  continuous  point  source  in  a  wind,  let  the  rate  of  emission  be 
g(t')  such  that  in  a  short  interval  from  t'  to  t'  +dt'  an  amount  qdt'  is 
emitted.  Each  of  these  ’puffs'  generates  its  own  cloud,  and  the  total 
concentration  field  is  obtained  by  the  summation  of  contributions  from 
the  individual  puffs.  For  a  source  maintained  indefinitely,  the  combined 
concentration  field  of  the  many  puffs  by  integration  is 


<3dt' 


(2n)  3/2  /DxDyDz  (t  -t’)3/2 


[x-U(t-t')  ]  2 _ _ z2 

2Dx(t-f)  2  Dy(t-t')  2  Dz{t-f) 


[4]  Crank,  J.,  The  Mathematics  of  Diffusion,  Oxford  University  Press, 
1956. 

[5]  Csandy,  G.  T.,  Turbulent  Diffusion  in  the  Environment,  Reidel  Pub¬ 
lishing  Go.,  Boston,  1972. 
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When  *2/Ox  >5>  5 r/D^  and  z2/Z>z 


the  solution  is: 


JSL 


2  n  UOyOfi 


exp 


f-  -4  -  -”i 

l  20J,2  2oz2J 


(2-5) 


Equation  (2-5)  is  the  plume  solution. 


The  dosage  for  the  puff  model  is  obtained  by  inserting  Equation 
(2-3)  into  Equation  (2-1)  and  integrating.  The  resulting  equation  is: 


(2-6) 

where  erf  is  the  error  function. 

The  plume  model  does  not  depend  on  time.  Hence,  the  dosage  in  a 
given  cell  is  simply  equal  to  the  product  of  the  nth  power  of  concentra¬ 
tion  and  the  total  evaporation  time  te  ?  that  is: 


v 


c 


(2-7) 


DERIVATION  OF  INTERIOR  DOSAGE 

The  time  variation  of  the  toxic  vapor  inside  a  building  is  propor¬ 
tional  to  the  infiltration  rate  and  the  difference  of  outside  and  inside 
vapor  concentration.  The  equation  is  written  as: 


dCj 

V  —  =  I(C0-Cr) 
dt  u  x 


(2-8) 


where  V  is  the  total  empty  space  of  the  building?  I  is  the  rate  of  in¬ 
filtration  (volume  per  unit  time);  Cq  is  the  outside  vapor  concentration; 
and  Cj  is  the  inside  vapor  concentration.  Equation  (2-8)  can  also  be 
written  as: 

+  RCj  =  RCq  '  (2-9) 
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where  R  *  I/V  is  the  specific  infiltration.  The  infiltration  depends 
mainly  on  the  tightness  of  the  construction  and  on  the  wind  velocity  [6J . 
For  high-rise  buildings,  the  chimney  effect  is  also  important.  In  the 
case  of  a  building  with  many  partitions  and  of  tight  construction,  air 
may  enter  on  the  windward  side  in  such  quantity  as  to  build  up  a  slight 
positive  pressure  and  thereby  reduce  infiltration.  In  general,  it  can 
be  assumed  that  the  air  which  enters  the  building  on  the  windward  side 
is  equal  to  that  which  leaves  on  the  leeward  side.  Empirical  equations 
for  infiltration  have  been  developed.  In  tne  VM,  the  following  equation 
is  employed  [2] : 


R  =  0.25  +  0.02165  V  0.00833  jAr| 


where  U  is  the  wind  velocity  in  miles  per  hour,  and  At  is  the  tempera¬ 
ture  difference  between  the  interior  and  exterior  of  the  building  in  °F. 

Multiplying  Equation  (2-9)  by  the  integration  factor 


and  integrating,  we  obtain 


t  t* 

-/  Rdt  t  J  Rdf 
=  e  0  RCQ  e  0  dt ' 

0 


Since  R  is  independent  of  time,  Equation  (2-10)  becomes: 


(2-10) 


„  -Rt  t  „  Rt'  , 
Ct  =  R  e  \  Cfi  e  dt 


(2-11) 


For  the  plume,  model,  the  outside  vapor  concentration  does  not 
depend  upon  time,  so  that 


Cj  ~  t? 


-Rt 


c 

f  Rf 

J  6 


dt'  = 


1SL 


2tt  t/o. 


f  -Rt 

1 

2  2 1 
f  z *■ 

r°  j 

exp 

"  2 

i  +  ~2 

°y 

(2-12) 


[6]  Jennings,  B.  H.,  Environmental  Engineering,  International  Textbook 
.  Company,  Scranton,  Pennsylvania,  1970. 
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for  which  Equation  (2-5),  the  continuous  spill  solution,  has  been  used 
for  Cq  .  Equation  (2-12)  shows  that  at  t  *0  ,  Cj  »0  and  that  at  t  =  *, 
C1  ~cO  •  The  vaP°r  concentration  inside  the  building  is  attenuated  by 
the  factor 


(l 


For  the  puff  model,  substituting  the  instantaneous  spill  solution, 
Equation  (2-3),  into  Equation  (2-11)  for  Cq  we  obtain: 


C  r  ~ 


2mR 


1  (2n)V20x0yo2 


cxn 


ai  +  £. 

°z 


-Rt 


<? 


c 


exp 


(x-Ut')2 


2oi 


+  Rt  * 


dt‘ 


2m  Roxf2v 


(2n)3/2  oxoyo2  2U 


exp 


1 

y2  z2 

xR  o2r  2' 

“  2 

2  +  T 

[°y  °lj 

exp 

.  U  2 t/2 

-Rt 


fer/f 


L 


^  ^  ax 


+  erf 


( Ut-lx+RaJ/U ) 


2  / »t\ 


/2  Ov 


(2-13) 


In  Equation  (2-1)  if  n-1  the  indoor  dosage  is: 


-f 

J  0 


vj  =  Cj  dt 
'0 


»  rt 

_  -Rt  .  Rt-'  .. 

Re  dt  Cq  (t)  e  dt 

i  Jo 


C0(t)dt  =  v0 


(2-14) 


If  n?<l  ,  then 


f  CO  #00 

Vj  =  Cj  dt  dt 

Jo  Jo 


=  Vn 


(2-15) 
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This  means  that  the  people  inside  and  outside  the  building  will  inhale 
the  same  amount  of  toxic  gas  but  will  suffer  different  effects  because 
of  the  time  delay  due  to  the  infiltration. 


The  general  equation  of  dosage  for  the  plume  model  is  obtained  by 
inserting  Equation  (2-12)  into  Equation  (2-1). 


VI  = 


1SL 


[2lU0yVZ J 


exp 


Oil  +  *1 

2  2 


*  n 


dt 


13L 


2tt  UOyOz 


exp 


(xl 

°v  +  °z 


t«-r 


n(n-l) (n-2)  (  -3Kte] 
3*3!  lA"e  > 


+  (— l)J^*~^>  —  (^**1)  (n— 2 ) . « .  (n-p-1 )  f  ^ 

p«p!  '•'L  J +  *  *  * 


(2-16) 


Here,  instead  of  using  the  upper  limit  of  integration  of  infinity,  a 
finite  number,  te,  the  evaporation  time,  is  used  and  the  binomial  series 
expansion  is  employed. 

Substituting  Equation  (2-13)  into  Equation  (2-1)  gives  the  indoor 
dosage  for  the  puff  model  as  follows: 


VI  = 


2m  /2tT  Rox~]n 

n 

f£  * *2  f 

xR  OvR2 

( 2ir)  3^2oxaya2  2  U  exP 

*\ 

tr  +  oT\ 
l- tj  az  ) 

exp 

-  +  -* — y 

U  2  V2  j 

foo 

e~nRt 

erf 

f  2  \ 

x+oxR/l7J 

—  erf 

'x-K)$R/U  -  yt] 

‘*0 

.  ^2*  ax  j 

dt 


(2-17) 
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Since  it  is,  in  general,  a  fractional  or  irrational  number,  no  closed 
form  of  solutions  of  Equation  ( 2 ~ 1 V )  is  expected.  Numerical  integra¬ 
tion  is  used  to  evaluate  Equation  ( 2 — 1 V )  .  'l'ho  computing  time  can  bo 
reduced  substantially  with  consideration  of  the  property  that  erf(p)  «1 
when  p  >  5  . 


DISCUSS  ION 

The  modification  of  the  plume  and  puff  models  for  toxic  gas  changes 
the  calculation  from  the  time-step  method  to  direct  integration.  It 
improves  the  accuracy  and  saves  computing  time  and  storage  space.  Since 
there  are  other  models  that  still  use  the  time-stop  method  for  calcula¬ 
tion,  it  is  impossible  to  change  VM  structure  at  this  stage  to  elim¬ 
inate  the  time-stop  procedure  entirely. 

Four  tost,  runs  have  been  made  to  compare  the  modified  VM  to  the 
original  VM.  Two  are  simulations  of  ammonia  spills  using  the  puff 
model,  and  the  other  two  are  simulations  of  chlorine  spills  using  the 
plume  model.  For  each  chemical,  one  of  tho  two  runs  is  with  the  orig¬ 
inal  VM  and  the  other  is  with  the  modified  VM.  The  following  table 
compares  tho  results  between  the  original  and  modified  runs. 


COMPARISON  OF  RFSULTS  URTMKN  ORIGINAL  AN!)  NODI  FIND  VN 


CHEMICAL 

ORIGINAL  VM 

MODIFIED  VM  I 

Indoors 

Outdoors 

Indoors 

Ammonia 

heaths 

5 , 560 

1,765 

5,298 

4,770 

Chlorine 

Deaths 

Injuries 

8 , 659 

3,670 

0,281 

0 

7,333 

5,531 

7,089 

6,390 

For  ammonia,  the  outdoor  toxic  deaths  for  both  runs  are  similar,  but 
tho  indoor  deaths  are  considerably  greater  for  the  modified  case.  No 
injuries  are  shown  because  no  probit  equation  exists  in  the  VM  for  toxic 
injuries  from  ammonia.  For  chlorine,  the  indoor  and  outdoor  toxic 
deaths  for  the  original  and  modified  cases  are  similar.  However,  for 
injuries,  the  modified  VM  yields  a  significantly  greater  number,  par¬ 
ticularly  indoors  whore  tho  old  VM  computed  no  injuries. 

Tho  input  data  for  these  test  runs  are  given  in  Appendix  A.  Figure 
A-l  is  a  map  of  the  geographical  area  showing  tho  spill  sir,e  and  loca¬ 
tion,  tho  wind  direction  and  speed,  and  tho  location  of  the  population 
colls.  Figure  A-2  presents  the  Goographic/Domographic  file  in  the  VM 
which  shows  the  population  and  vulnerable  characteristics  of  each  popu¬ 
lation  cell.  (Note  that  coll  numbers  R31,  R32,  etc.,  are  river  cells 
which  are  not  shown  on  the  map.)  Figures  A- 3  and  A-4  give  tho  input 
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data  Cor  the  ammonia  and  chlorine  spills,  respectively,  and  Figure  A-5 
(extracted  from  reference  (7] ^presents  the  dictionary  for  interpreting 
the  input  field  numbers.  The  modified  program  saves  about  20%  in  com¬ 
puting  expense  including  computing  time  and  storage  over  the  original 
program. 
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NOMENCLATURE  AND  UNITS 


C  vapor  concentration 

CQ  outdoor  vapor  concentration 

Cj  indoor  vapor  concentration 

Dx,Dy,Dz  diffusivity 

I  infiltration 

q  spill  rate 

m  total  spill  mass 

R  specific  infiltration 

t  time 

U  wind  velocity 

v  dose 

V  volume 

x,y,x  coordinates 

ox,Oy,o2  dispersion  coefficients 


g/cm3 

g/cm3 

g/cm3 

cmi/s 

cm3/s 

g/s 

9 

s“l 

s 

cm/s 

cm3 

cm 

cm 


(7J  Rowley,  R.  M.,  and  A. 
Enviro  Control,  Inc., 


H.  Rausch,  Vulnerability  Model  User's  Guide, 
October  1977. 


11 


Chapter  3 

MODIFICATION  OF  THERMAL  INJURY 
AND  LETHALITY  CRITERIA  FOR  HUMANS 

INTRODUCTION 

In  the  Vulnerability  Model,  the  fire  damage  to  personnel  is  assessed 
by  the  probit  equation  which  is  defined  as  [lj: 

Pr  =  a  +  b  in  v 

where  v  is  the  thermal  dosage  and  a  and  b  are  determined  from  existing 
experimental  data.  The  thermal  dosage,  v,  depends  on  the  thermal  radi¬ 
ation  intensity,  I,  and  the  exposure  time,  ft.  In  general,  the  thermal 
dosage  can  be  written  as: 

v  =  tla 

where  the  index  o  is  also  determined  from  experimental  data.  In  the  VM, 
coefficients  a  and  b  and  index  a  had  been  based  on  nuclear  explosion 
data  [8,9].  The  thermal  radiation  from  nuclear  explosions  is  primarily 
in  the  ultraviolet/  visible  spectra  with  high  intensity  and  short  dura¬ 
tion.  On  the  other  hand,  the  thermal  radiation  from  hydrocarbon  com¬ 
bustion  is  in  the  infrared  range  with  low  intensity  and  longer  duration. 
Butterfield  [10]  has  shown  that  almost  twice  as  much  visible  radiation 
is  required  to  produce  an  equivalent  injury  caused  by  an  infrared  radi¬ 
ation  source.  Owing  to  these  differences,  data  from  sources  other  than 
nuclear  explosions  have  now  been  used  to  determine  the  coefficients  a 
and  b,  and  the  index  a. 


[8]  Glasstone,  S.  (ed.),  The  Effects  of  Nuclear  Weapons,  USAEC,  April 
1962. 

[9]  White,  Clayton,  S.,  The  Nature  of  the  Problems  Involved  in  Estimat¬ 
ing  the  Immediate  Casualties  from  Nuclear  Explosions ,  CEX  71.1, 
Lovelace  Foundation  of  Medical  Education  and  Research,  Albuquerque, 
New  Mexico,  March  1971. 

[10]  Butterfield,  W.  J.  H.,  E.  R.  Drake  Seager,  et  al.,  Flash  burn  from 
atomic  weapons,  Surgery,  Gynecology  and  Obstetrics  103 (6) : 655-665, 
December  1956. 
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DERIVATIONS 
Exponent  a 


*®>^=eg p 


In  the  VM,  the  thermal  dosage  is  expressed  by 

v  «  tl1*/3 


(3-1) 


where  t  is  in  seconds  and  I  is  in  joule/m?-sec.  The  data  from  which  the 
value  of  c.  -4/3  was  deduced  are  reproduced  in  Table  3-1.  It  is  obvious 
that  Equation  (3-1)  is  for  fire  deaths.  However,  the  same  expression 
had  been  used  in  the  VM  for  fire  injuries. 


Table  3-1 .  RELATIONSHIP  OF  DEATH  FROM  RADIATION  BURNS 
TO  RADIATION  LEVEL  AND  DURATION 


Portion 

Killed 

<%> 

Duration 

(s) 

Radiation  Intensity 

Dosaqe 

tiV3 

(cal/cm2-s) 

( joule/m2 -s) 

1 

1.43 

3.50 

146,000 

1099  x  xo4 

1 

10.1 

0.792 

33,100 

1073  x  IQ4 

1 

45.2 

0.243 

*  10,200 

1000  x  104 

50 

1.43 

6.30 

263,600 

2417  x  lo4 

50 

10.1 

1.385 

57,950 

2264  x  lo4 

50 

45.2 

0.442 

18 , 500 

2210  x  104 

99 

1.43 

14.0 

586,000 

7008  x  104 

99 

10.1 

3.07 

128,000 

6546  x  104 

99 

45.2 

0.952 

39,800 

6149  x  lo4 

The  sequence  of  data  is  for  20-KT,  1-MT,  and  20-MT  weapons. 
The  dosage  has  been  calculated  for  the  radiation  intensity 
in  joule/m2 -s  , 


Figure  3-1  shows  the  human  skin  tolerance  timo  to  absorbed  thermal 
energy.  The  blister-line  corresponds  to  second-degree  burns  and  the 
survival-line,  to  first-degree  burns.  The  slope  of  the  pain-lino  is 
-4/3,  that  of  the  survival-line  is  -1.35,  and  that  of  the  blister-line 
is  -1.375.  Therefore,  the  thermal  dosage  for  pain  is 

v  »  tl4/3 

and  for  first-degree  burns  is 

Vj»  tJ1*35 

and  for  second-degree  burns  is 

Vj“  tl1,375 
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'RANCE  TIME  (Ste.J 


There  are  no  curves  given  for  third-degree  burns  and  fire  deaths.  In 
view  of  these  values,  a  =  4/3  has  been  retained  for  the  thermal  lethal¬ 
ity  and  injury  (first-degree  burns)  probit  equations  in  the  VM. 


g  it  innnB 
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ABSORBED  ENERGY  RATE  (CAL 


Figure  3-1.  Human  Skin  Tolerance  Time  to  Absorbed  Thermal  Energy 
Delivered  in  a  Rectangular  Heat  Pulse 


Thermal  Injury  Probit 

In  the  VM,  the  probit  equation  for  thermal  injury  was: 

Pr  =  -42.25  +  3.0186  In  (t I4/3 )  (3-2) 

When  injury  is  1%  (Pr  =  2.67),  the  dosage  is: 

v  =  tI4/3=  2,902,500 

In  Figure  3-1,  the  dosage  for  first-degree  burns  is  calculated  as: 

tl4/3  =  1,280,000 

In  reference  (10) ,the  threshold  for  first-degree  burns  is  0.94  cal/cm2 
(3.93  *  104  joule/m2)  for  a  one-second  duration.  Thus: 

tl4/3  =  1,. 336, 100 

In  reference  [11], the  threshold  for  first-degree  burns  is  1.2  cal/cm2 
(5.021  x  104  joule/m2)  for  a  duration  of  3  seconds.  Thus: 

tl4/3  =  1,284,300 

The  average  of  the  last  three  figures  is  approximately  1,300,000.  The 
value  of  2,902,500  from  the  VM  is  about  2.23  times  higher  than  this 
average.  Using  the  factor  of  2.23  to  modify  Equation  (3-2),  we  obtain 
a  new  probit  equation  for  thermal  injury: 

Pr  =  -39.83  +  3.0186  In  (tl4/3)  (3-3) 

Although  this  adjustment  factor  has  been  derived  on  the  basis  of  the  1% 
level  of  first-degree  burns,  it  is  assumed  to  hold  for  all  levels,  in 
the  absence  of  specific  data  for  these  other  levels. 

Thermal  Lethality  Probit 

The  probit  equation  for  burn  deaths  in  the  VM  was: 

Pr  =  -14.9  +  2.56  In  (tl4/3/ 104) 

For  burn  deaths,  no  other  data  besides  the  nuclear  explosion  are  avail¬ 
able.  Using  the  same  modifying  factor  of  2.23,  the  probit  equation  for 
deaths  now  becomes: 

Pr  =  -12.8  +  2.56  in  (tJ4/3/104) 


[11]  Hardy,  J.  D.,  H.  G.  Wolff,  and  H.  Goodell,  Studies  in  pain,  new 
method  for  measuring  pain  threshold:  Observations  on  spatial  sum¬ 
mation  of  pain.  Journal  of  Clinical  Investigation  19,  1940. 
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NOMENCLATURE  AND  UNITS 
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radiation  intensity 


probit 
thermal  dose 
index 


J/m  -s 
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Chapter  4 

REMODELING  OF  FLASH  FIRE 


INTRODUCTION 


The  flash  fire  .model  had  been  designed  using  assumptions  and  data 
deemed  appropriate  for  the  first  stage  of  development.  Some  of  the 
simulated  results  were  incompatible  with  experimentally  observed  data. 
Consequently,  the  model  may  underestimate  injuries  and  deaths  under  cer¬ 
tain  conditions.  There  are  three  assumptions  made  during  this  first 
stage  of  development  that  should  be  modified  to  correct  for  these  defi- 
ciences.  First,  in  the  flash  fire  model  the  air-fuel  mixture  is  divided 
into  two  parts  [1] .  The  concentration  of  the  outer  part  varies  from  the 
lower  flammable  limit  to  stoichiometric.  The  concentration  of  the  inner 
part  varies  from  stoichiometric  to  upper  flammable  limit.  It  is  assumed 
that  the  outer  portion  burns  completely  and  that  the  inner  portion  burns 
incompletely  from  lack  of  sufficient  oxygen  content.  This  assumption 
underestimates  the  vapor  quantity  that  takes  part  in  combustion  and, 
in  turn,  underestimates  the  total  thermal  radiation.  For  example,  in 
Appendix  C2  of  the  VM  [1] ,  the  ratio  of  upper  flammable  limit  to  lower 
flammable  limit  of  methane  is  2.642  ( [1] ,  Table  C2-1) .  From  Figure 
C2-1  in  reference  II] ,  it  is  found  that  the  amount  of  fuel  consumed  is 
about  40%  of  the  total  vapor.  The  real  situation  is  that,  after  the 
fire  starts,  the  hot  gas  causes  a  strong  turbulent  flow  to  speed  up  the 
mixing  of  the  fuel  with  the  ambient  air.  A  small  percent  of  the  fuel 
vapor  will  escape  to  the  ambient  atmosphere,  but  the  remainder  will  burn 
out. 

The  second  assumption  is  that  the  combustion  is  instantaneous. 
Therefore,  the  thermal  radiation  from  the  burning  fuel  is  neglected. 
Actually,  the  thermal  radiation  in  the  combustion  period  is  of  the  same 
order  as  that  in  the  cooling  period.  The  total  combustion  time  of  the 
fuel  vapor  depends  on  the  mixing  process  and  the  flame  speed. 

Third,  in  the  cooling  period,  the  hot  gas  is  assumed  to  be  cooled 
from  the  initial  temperature  T±  to  a  "half  life"  temperature  Tg ,  which 
is  defined  as  the  average  of  the  initial  temperature  Ti  and  the  ambient 
temperature  Tq,  i.e.,  Tg  =  (T±  +  Tq)/2  .  Based  on  this  assumption,  an 
expression  for  cooling  time  t  is  derived  (Equation  (4-18),  reference  II]) 
The  cooling  time  computed  from  this  equation  is  improper.  For  example, 
the  cooling  time  for  a  gas  with  a  higher  initial  temperature  or  a  larger 
vapor  mass  is  shorter  than  that  for  a  gas  with  a  lower  initial  tempera¬ 
ture  or  a  smaller  vapor  mass.  The  "half  life"  criterion  which  is  used 
in  linear  problems  such  as  nuclear  radiation  decay  and  mechanical  vibra¬ 
tion  damping  is  not  suitable  for  a  nonlinear  problem  such  as  thermal 
radiation. 


To  remedy  this  difficulty  and  to  base  the  flash  fire  simulation  more 
accurately  on  physical  reality,  a  new  flash  fire  model  has  been  developed 
using  a  modification  of  the  fireball  model  [12,13,14). 


ASSUMPTIONS  IN  THE  MODIFIED  FIREBALL  MODEL 

The  fireball  model  was  originally  designed  for  liquid  propellant 
rocket  explosions  on  a  launchpad.  The  general  assumptions  for  that 
model  are  as  follows  (13,14): 

(1)  The  rate  of  fuel  addition  to  the  fireball  is  constant. 

(2)  A  stoichiometric  mixture  is  assumed  to  exist  at  ignition. 

(3)  All  the  available  fuel  participates  in  the  reaction. 

(4)  The  fireball  is  an  isothermal,  homogeneous  body  which  is 
spherical  at  all  times. 

(5)  The  fireball  radiates  as  a  blackbody. 

(6)  The  fuel  burnout  time  and  the  fireball  liftoff  time  coin¬ 
cide. 

Assumptions  (1)  and  (6)  can  only  be  applied  to  a  liquid  fuel  spilled  on 
the  ground  and  are  inapplicable  to  the  present  problem.  Assumptions  (2) 
to  (5)  are  applicable  to  the  modified  fireball  model.  One  additional 
assumption  for  the  present  fireball  model  is  that  the  flame  speed  is 
constant,  A  laboratory  experiment  [15)  showed  that,  as  a  first  ap¬ 
proximation,  the  flame  propagation  took  place  at  nearly  constant  veloc¬ 
ity.  Another  assumption  which  was  used  but  not  mentioned  in  the  origi¬ 
nal  fireball  model  is  that  the  fire  started  at  the  center  and  propagated 
radially  outward.  We  will  retain  this  assumption  also. 

The  general  assumptions  in  the  modified  flash  fire  model  are: 


(12)  Kite,  F.  D.,  and  B.  E.  Bader,  Pad-Abort  Thermal  Flux  Model  for 
Liquid  Rocket  Propellants,  SC-RR-66-577 ,  Sandia  Laboratory,  Albu¬ 
querque,  New  Mexico,  November  1966. 

(13)  Bader,  B.  E.,  A.  B.  Donaldson,  and  H.  C.  Hardee,  Liquid-propellant 
rocket  abort  fire  model.  Journal  of  Spacecraft  and  Rockets  8(11): 
1216-1219,  1971. 

(14)  Hardee,  H.  C.,  and  D.  0.  Lee,  Thermal  hazard  from  propane  fireballs, 
Transportation  Planning  and  Technology  2: 121-128,  1973. 

(15)  Leyer,  J.  C.,  C.  Guerraud,  and  N.  Manson,  Flame  propagation  in 
small  spheres  of  unconfined  and  slightly  confined  flammable  mix¬ 
tures,  in  Fifteenth  Symposium  on  Combustion,  The  Combustion  Insti¬ 
tute,  Penn  State  University,  1974. 
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(1)  Fire  starts  at  the  center  and  propagates  with  constant 
velocity  outward. 

(2)  The  fuel-air  mixture  is  stoichiometric  and  homogeneous. 

(3)  The  fireball  is  an  isothermal  spherical  body. 

(4)  The  fireball  is  a  blackbody. 

(5)  All  the  fuel  participates  in  the  reaction. 


DERIVATIONS 


Burning  Mass 


Consider  a  fireball  of  radius  r{t).  The  flame  front  propagates  out¬ 
ward  and  advances  a  distance  $r  in  the  time  St.  The  energy  balance 
equation  in  the  thin  layer  between  r  and  r  +  Sr  is: 


p£tf4irr25r  =*  T.  PiCpj  (T  -To)47rr25r  +  eo(T4-  TQ)4irr26t 


(4-1) 


where:  pf  is  the  fuel  density;  H  is  the  heat  of  combustion;  is  the 
density  and  Cpi  is  the  specific  heat  at  constant  pressure  of  the  ith 
product;  Tg  is  the  ambient  temperature;  a  is  the  Stefan-Boltzmann  con¬ 
stant;  and  e  is  the  emissivity  of  the  gas  which  is  equal  to  unity  in  the 
present  case.  Let  Sr/St  *  S  ,  the  flame  velocity.  Then  Equation  (4-1) 
becomes : 


spfH  -  s  £PjCpi  (r-r0)  +  eo(Tk- r04) 


(4-2) 


The  combustion  is  isobar ic»  so  that  the  density  p  depends  upon  tempera¬ 
ture  T  only ,  or  p  *  p(T)  ,  The  specific  heat  Cp  is  also  a  function  of 


temperature i 


Therefore,  when  the  flame  velocity  5  is 


known,  the  flame  temperature  T  can  be  solved  from  Equation  (4-2) .  On 
the  other  hand,  if  the  flame  temperature  is  measured,  then  the  flame 
velocity  is  obtained  from  Equation  (4-2).  Usually,  in  laboratory  ex¬ 
periments  the  flame  velocity  is  measured,  and  in  field  tests  the  flame 
temperature  is  measured  (16,17J . 

Due  to  turbulent  mixing,  most  of  the  gas  in  a  vapor  cloud  will  be 
burned  once  ignition  takes  place.  As  a  conservative  assumption,  we  will 
assume  that  the  portion  of  the  vapor  cloud  with  concentration  above  the 
lower  limit  of  inflammability,  C^,  will  be  burned  in  the  flash  fire. 

(16]  Strehlow,  R.  A.,  Fundamentals  of  Combustion,  International  Textbook 
Co.,  Scranton,  Pennsylvania,  1968. 

(17)  American  Gas  Association,  LNG  Safety  Program.  Interim  Report  on 
Phase  II  Work,  July  1974. 
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By  rewriting  Equation  (2-3)  as 


2  m 

C  “  (2W)  V^0x0y02  SX^ 


f  (x-Pt 
1“  20J 


202  2 o|  j 


and  letting  (x-l/t)/ox  =  x*  ,  y/oy  =  y’  ,  z/o2  =  z'  ,  we  obtain: 


(4-3) 


C  = 


2  m 


(2ii)  V2oxoyoz 


exp 


[- i  (,■»*„•’♦  «->y 


_  2m 

(2  n)  i/2axOyOz 


exp 


•H 


(4-4) 


Equation  (4-4)  transforms  the  ellipsoid  into  a  sphere  with  r‘  - 
2  +  y'2  +  z'2.  When  C  =  CL  ,  from  Equation  (4-4) 


rL  " 


2  in 


2m 


(2n)  Vs 


GX°UC2CL 


r 

j. 


»/2 


(4-5) 


is  the  radius  of  equi-concentration  surface.  The  total  fuel  consumed 
in  the  fire  is  the  mass  inside  the  sphere  of  radius  . 


m 


Mf  (2n)  3/2 


Px°y°z 


exp 


dx  dydz 


m 

~  (211)  3/2 


x,2  +  y,2  +  Z,z 
2 


dx'dy'dz' 


(4-6) 


To  carry  out  the  integration,  it  is  better  to  transform  Equation  (4-6) 
from  Cartesian  coordinates  (x',y',2')  to  spherical  coordinates  (r,6.<f>) 
where  0  is  the  latitude  and  <f>  is  the  azimuthal  angle.  Under  this  trans¬ 
form,  Equation  (4-6)  becomes 


Mf  = 


ro 

(2 it)  3/2" 


rL 
■  0 


\  _2 


-;r 


r2dr 


n 

sin  6diJ> 
0 


211 

d4> 

jO 


(4-7) 
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For  a  real  spill,  erf (r^Z/2  )  -  1  and  rL  e  rL  2  =*  0,  then  Mf=m,  the  total 
mass  of  fuel.  For  the  plume  model,  using  Equation  (2-5)  and  following 
the  same  procedure,  it  will  turn  out  that  Mf  *  qte  ,  where  te  is  the 
total  evaporation  time.  The  maximum  fireball  radius  r £  is  determined 
from  the  total  mass  of  products,  Mq  =  Mf  +  ,  where  Ma  is  the  mass  of 

air, 

4 


Nq  =  2  Op 


(4-8) 


where  p p  is  the  average  density  of  products  at  temperature  Tf .  For 
ppTf  *  Pporo  where  Pp0  is  the  average  density  at  ambient  temperature 
T0,  thens 


= 


3  Tf  )  >/J 


4  71 


Ppo  To 


(4-9) 


Since  the  fuel  mass  Mf  is  given,  it  is  easier  to  estimate  the  radius  r*> 
by  the  following  relation, 


rb 


_  r 3  Tf  i 

[4  7Tpf0  To  J 


1/3 


(4-9a) 


where  pfo  is  the  fuel  density  at  ambient  temperature. 
Burning  Phase  Thermal  Dose 

The  probit  equation  for  thermal  damage  to  personnel  is; 

Pr  *  a  +  b  In  v 


(4-10) 


where  a  and  b  are  constants  and  v  is  the  thermal  dosage.  The  thermal 
dosage  depends  upon  the  radiation  intensity  absorbed  by  the  body  and 
the  duration.  In  the  VM,  the  thermal  dosage  is  given  by  the  equation: 

v  =  (4-11) 

The  radiation  intensity  is  calculated  from  the  equation  as  follows: 


I  =  acoFuT 


(4-12) 


where  a  is  the  absorptivity  of  the  body  and  F\z  is  the  view  factor. 
The  view  factor  for  a  flat  surface  of  unit  area  at  a  distance  d  from 
the  fireball  center  is  given  by  (18} : 


F\z 


(4-13) 


(18)  Love,  T.  J.,  Radiative  Heat  Transfer,  Merrill  Publishing  Company, 
Columbus,  Ohio,  1968. 
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With  Equations  (4-12)  and  (4-13),  the  thermal  dosage  becomes: 


-ft 


2 


(4-14) 


For  a  growing  fireball,  the  flame  temperature  T  -Tf  is  constant  and 
dt  -dr/3  ,  Equation  (3-14)  can  be  integrated  as: 

1  t  h\ <4/3  3  1  f  ^4/3  TM8/3 

V  -  -  (oeorf)  Jo  (?)  dr  =  -  -  (neo Tf  j  '  [~J  rb  (4-15) 


Cooling  Phase  Thermal  Dose 

After  the  fuel  is  burned  out,  the  fireball  starts  to  cool  down 
through  the  radiation  heat  loss.  The  energy  equation  of  the  cooling 
process  is  in  the  form: 


-  I  Cp^  V  dT  =  A£o(t4-  Tg)dt 


(4-16) 


where  V  is  the  volume  and  A  is  the  surface  area  of  the  fireball. 
Equation  (4-16)  we  have: 


‘1  picPi 


d t  =  -  4 — -77  dT 

A£0(T4-  To  ) 


(4-17) 


and  by  substituting  Equation  (4-17)  into  Equation  (4-14),  we  obtain: 


Te  Picpi 

p.-taW"  \  v 

dzj  A (T  -  T»  ) 


(4-18) 


where  the  integration  limit  Te  is  the  final  temperature.  Since  V/A  =  r/ 3, 
r  =zi>  (T/Tf ) l/3  and  pjT  =  PioTo>  Equation  (4-18)  can  be  transformed  to: 


V  =  (a4£0)  l/3  T„  ^ 


,rhl»/«_£fe_  r“/,?pi.ppi  . 

=  H  3Tf'v»  _  T-.T>  " 


(4-19) 


The  specific  heat  Cpi  is  a  function  of  temperature.  Empirical  equations 
are  designed  to  fit  the  experimental  data.  However,  the  best-fit 


equations  are  not  always  in  a  similar  form  [19].  For  example,  the  equa¬ 
tions  for  air,  hydrogen,  carbon  dioxide,  and  steam  are  as  follows: 


air: 

H2 : 
C02: 
H20: 


Cp  »  a i  +  b\T  +  C\T 
Cp  =  a2  +  b2T  +  c2/Tl/2 
Cp  =  a 3  +  bi/T  +  c3/T2 
Cp  =  ah  +  2?4/rl/4+  c4/T 


where  the  a's,  b's,  and  c's  are  constants.  For  a  first  approximation, 
we  will  assume  that  Cp  is  constant.  Integration  of  Equation  (4-19) 
yields: 


v  =  9(a4eo) l/3  (£ 
i 

Pi  o  cpi  )  To  |  ^ 

I8/3  rb 

1  rp9" 

•  { = 

1 

m  23  /9  ?  1 

V 

-13  -36 (j-1) 

[TfJ  e  fit  “13  -36 (j-1) 

[reJ 

- 

(4-20) 

Here  the  series  expansion  of  [1  -  (To/T)k )~l  has  been  used.  In  general, 
the  final  temperature  Te  is  the  ambient  temperature  T0  .  Here  we  will 
choose  Te  =  1.1  Tq .  This  choice  will  not  affect  the  results  of  thermal 
damage,  because  at  such  a  low  temperature  the  radiation  effect  is  neg¬ 
ligible,  but  it  will  enhance  the  convergence  of  the  series.  When 
T0/Te  =  1/1.1  ,  the  series  is  equal  to: 


r  =  - ± _ 

j=l  -13  -36  (j-1) 


1  W 


1.1 


=  -.0263656 


and  then: 


v  -  9(a**ea) 1/3  (Epi0cpi )r0 L*’!  n/9 

[ t|3/9  E  -rrr — :  T r  f~]  J+  .02067  T( 
f  j= o  -13  -36(j-l)  (Tfj  1 


(4-21) 


[19]  Faires,  V.  M.,  Thermodynamics,  5th  edit.,  Macmillan  &  Company,  New 
York,  1970. 


23 


In  the  VM,  tho  assessment  of  fire  damage  to  structures  is  based  on 
studies  of  the  ignition  of  wood  [1] .  Factors  influencing  wood  ignition 
ares  (1)  radiation  intensity  level,  (2)  duration  of  radiation  exposure, 
and  (3)  wood  type.  Wood  typo  is  not  treated  explicitly;  average  values 
are  used.  The  criteria  for  ignition  of  structures  are; 

1.  The  radiation  intensity,  An  '  must  exceed  the  value: 


I„,  =  1.34  x  io' 


erg/cm' -s 


(4-22) 


The  effective  duration  of  tho  radiation,  t eff>  must  exceed 
the  value  given  by: 


'■m 


f  1.23  v  10fl  erg/cm* -s ]  v;  1.364  x  1QU' 

l  I  "An  )  =  U  -I,,,)  V? 


(4-23) 


Tho  criterion  (4-22)  can  be  used  in  the  present  case.  However,  the 
criterion  (4-23)  has  to  be  modified  because  tho  radiation  intensity  I 
is  a  variable.  Wo  will  define  a  thermal  dosage  for  structural  ignition 

itS  J 


v  *  Ju  -rm) )/2dt 

If  v  >  1.364  *  10 1? ,  then  there  is  ignition  in  structures. 


(4-24) 


Figure  4-1  is  a  sketch  of  radiation  intensity  I  vs.  time  t.  To 
study  the  structural  ignition,  we  have  to  locate  the  points  l  and  2  and 
then  integrate  Equation  (4-24)  from  tt  to  t-u  and'  from  tj.  to  t2.  First, 
for  the  burning  fireball: 


An 


=  aco 


where  rm  is  tho  fireball  radius  corresponding  to  the  radiation  intensity 
An .  Substituting  An  into  Equation  (4-24),  we  obtain: 


v\ 
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vb 

dr 
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|  3/2 
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4S 

d2  j 

M  mJ  2  1 

+  l£m  in  rb  +(rl) 
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4) 


1/2 


(4-26) 
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Figure  4-1 .  Variation  of  Radiation  intensity  of  Fireball 


For  a  cooling  fireball,  the  thermal  radiation  intensity  is: 


r)2^ _ frhl 2  t14/* 


I  »  aeoKj  T  =  aeo 


~fpT 


Here  the  relation  r/rb  -  (T/Tf)  has  been  used. 
When  I  •  Im,  the  temperature  at  point  2  is: 


t,  .2-2/3) Vi- 


The  dosage  v2  is  obtained  as  follows: 
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(4-27) 


In  deriving  this  equation,  the  binomial  expansion  and  Equations  (4-17) 
and  (4-26)  are  used,  v  =  v\  +  v-i  is  the  total  radiation  dose.  If  v  is 
greater  than  or  equal  to  vm  ( vm  =  1.364  x  lO**),  the  structure  will  be 
ignited . 


DISCUSSION 

The  modification  of  flash  fire  is  based  on  the  fireball  model.  For 
both  the  puff  and  the  plume,  it  has  been  shown  that  the  fuel  consumed  in 
the  fire  is  very  close  to  the  total  evaporated  mass.  The  equivalent 
fireball  diameter  is  determined  from  the  total  evaporated  mass  (Equation 
(4-8)).  For  the  case  of  vapor  clouds  arising  from  instantaneous  spills, 
the  spherical  fireball  model  is  a  good  approximation  to  the  ellipsoidal 
puff  cloud.  However,  for  continous  spills  which  form  a  plume,  the  shape 
of  a  plume  is  an  elongated  truncated  ellipsoid.  When  the  ignition  cell 
is  not  far  from  the  spill  center,  the  spherical  fireball  is  a  reasonable 
approximation  for  the  plume  fire,  as  has  been  assumed  in  references  (13) 
and  (14).  But  when  the  ignition  cell  is  far  downwind,  the  fireball  model 
is  not  suitable  for  the  plume  because  of  the  significant  geometrical 
difference.  This  difference  affects  both  the  area  covered  by  the  vapor 
cloud  and  the  view  factor  calculation.  In  this  case,  either  an  ellip¬ 
soidal  or  a  cylindrical  fireball  model  should  be  used.  Neither  of  these 
has  been  developed  at  this  time.  Consequently,  for  a  continuous  spill, 
the  flash  fire  computation  is  presently  limited  to  situations  in  which 
ignition  occurs  fairly  close  to  the  spill  site  and  the  ratio  of  the 
major  to  minor  axis  for  the  truncated  ellipsoid  is  less  than  three. 
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In  the  energy  balance  equation  (4-1),  the  conduction  and  convective 
heat  transfer  are  not  included.  The  conduction  heat  transfer  is  smaller 
than  the  radiation  heat  transfer.  In  reference  [13J  the  authors  mention 
that,  for  a  flat  plate,  the  convective  heat  transfer  is  about  1/32  of 
the  radiation  heat  transfer.  For  a  fireball,  there  is  no  solid  boundary, 
and  it  is  by  no  means  clear  that  the  convective  heat  transfer  can  be 
neglected.  Consequently,  further  investigation  of  the  role  of  convec¬ 
tive  heat  transfer  in  the  fireball  model  is  needed. 

In  the  derivation  of  Equation  (4-1),  it  is  assumed  that  the  fireball 
is  an  isothermal,  homogeneous  body.  In  reality,  the  temperature  at  the 
fireball  center  is  the  highest  and  then  decreases  toward  the  boundary. 

The  temperature  obtained  from  Equation  (4-1)  for  a  given  flame  speed  is 
an  average  value  between  the  maximum  and  the  surface  temperature.  But 
in  the  damage  assessment,  the  thermal  radiation  is  from  the  fireball  sur¬ 
face.  In  Equation  (4-15)  the  thermal  dose  depends  upon  the  flame  speed 
and  fireball  surface  temperature.  Therefore  if  the  average  temperature 
obtained  from  Equation  (4-1)  is  used  in  the  dose  calculation,  an  over¬ 
estimated  damage  will  result. 

Take  methane  as  an  example.  The  following  table  is  obtained  from 
Equation  (4-1) 


Flame  temperature,  °K 

2000 

1900 

1800 

1400 

Flame  speed,  m/s 

4.44 

3.12 

2.20 

0.74 

From  laboratory  and  field  measurements,  the  flame  speed  is  in  the  range 
of  3  to  5  m/s  [20,21] .  The  field  measurements  of  LNG  fires  indicate  that 
the  surface  temperature  lies  between  1000  and  1400°K  [22] .  If  we  use  a 
value  of  4.44  m/s  for  the  flame  speed,  then  we  obtain  from  Equation  (4-1) 
a  flame  temperature  of  2000°K  and,  in  turn,  a  very  high  thermal  dose  from 
Equation  (4-15).  Similarly,  if  we  use  a  value  of  1400°K  for  flame  tem¬ 
perature,  then  we  obtain  a  lower  than  usual  flame  speed  of  0.74  m/s  from 
Equation  (4-1)  and,  in  turn,  from  Equation  (4-15)  we  obtain  a  very  high 
thermal  dose.  No  matter  which  measured  value,  temperature  or  flame 
speed,  is  used  in  Equation  (4-1),  an  unusually  high  thermal  damage  will 
result.  To  resolve  this  difficulty  at  the  present  time,  we  will  not  use 
Equation  (4-1);  instead,  we  will  use  nu  isured  values  for  both  flame  tem¬ 
perature  and  flame  speed  in  Equation  (4-15) . 


[20]  Lewis,  B.,  and  Von  Elbe,  Combustion,  Flames  and  Explosions  of  Gases, 
Academic  Press,  1951. 

[21]  Pangor,  E.,  Flame  Photometry ,  Van  Nostrand,  1967. 

[22}  Kanury,  A.  Murty,  Introduction  to  Combustion  Phenomena,  Golden  and 
Breach,  1975. 
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The  original  fireball  model  was  designed  for  rocket  fuel  and  oxi¬ 
dizer.  Since  there  is  no  nitrogen  in  the  fuel-oxidizer  mixture,  the 
flame  speed  is  higher  than  in  a  fuel-air  mixture.  Because  of  this  dif¬ 
ference,  a  comparison  of  the  results  of  the  original  fireball  model  with 
those  of  the  modified  flash  fire  model  would  be  meaningless. 

The  modified  flash  fire  model  is  divided  into  two  phases:  the 
burning  phase  and  the  cooling  phase.  The  duration  of  each  phase  depends 
on  the  flame  velocity  and  the  vapor  quantity.  The  burning  time  may  take 
several  minutes  for  large  vapor  clouds,  and  the  cooling  time  much  longer. 
In  an  actual  situation,  as  soon  as  the  fire  starts  the  exposed  people 
will  attempt  to  evade  the  radiation  by  running  to  a  shielded  location. 
Those  who  survive  or  are  not  seriously  injured  during  the  time  it  takes 
to  reach  a  sheltered  position  will  most  likely  not  be  killed  or  injured 
during  the  remainder  of  the  radiant  period.  Because  of  this  evasive 
action,  the  casualty  assessment  for  flash  fire  should  be  made  for  a 
radiation  intensity  integration  time  corresponding  to  the  time  required 
for  the  people  to  reach  shielded  positions.  In  a  normal  residential  or 
urban  area,  30  seconds  would  be  a  reasonable  time  for  most  people  to 
attain  shelter.  However,  for  people  located  on  a  beach  or  in  a  stadium, 
a  much  longer  integration  time  would  be  required.  Consequently,  the 
user  must  select  the  integration  time  to  suit  the  particular  situation 
being  studied.  For  the  case  of  structural  ignition,  the  integration 
time  would  be  for  the  entire  time  period  that  the  radiation  intensity 
exceeds  the  threshold  value. 

Four  test  runs  have  been  made  with  the  modified  flash  fire  model. 
These  are  made  for  an  instantaneous  LNG  spill  of  25,000  m3  at  the  same 
location  and  for  the  same  wind  conditions  as  shown  in  Figure  A-l 
(Appendix  A).  Thus,  the  geographic/demographic  file  of  Figure  A-2 
applies.  The  four  runs  correpond  to  integration  times  of  30  seconds, 

60  seconds,  entire  burning  time  (109  seconds),  and  infinity  (entire 
burning  and  cooling  time) .  The  table  below  presents  the  casualties 
for  these  four  cases.  Figure  A-6  in  Appendix  A  is  a  listing  of  the 
inputs  for  these  test  runs. 


Integration  Time 

Injuries 

Deaths 

Buildings  Destroyed 

30  s 

55 

30 

NA 

60  s 

1,969 

1,554 

NA 

109  s  (burning  time) 

5,808 

5,622 

NA 

Infinity 

10,374 

47,362 

5,648 

These  results  show  the  strong  dependency  of  tne  casualties  on  the 
integration  times.  Thus,  it  is  important  in  assessing  casualties  from 
flash  fire  to  use  realistic  estimates  of  time  required  for  evasive 
action. 
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NOMENCLATURE  AND  UNITS 


A  area 

Cp  specific  heat 

d  distance 

H  heat  of  combustion 

I  radiation  intensity 

n  index 

v  dosage 

r  radius 

S  flame  velocity 

t  time 

T  temperature 

V  volume 


cm2 

erg/g 

cm 

erg/g 
erg/cm2 -s 


cm/s 

s 


°K 


cm 


3 


«  absorptivity 

K  emissivity 

p  density  g/cm^ 

o  Stefan-Boltzmann  constant  erg/cm2-s-°K4 
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Chapter  5 

SUMMARY  OF  REVISIONS  TO  VULNERABILITY  MODEL  COMPUTER  PROGRAM 


INTRODUCTION 

Major  programming  revisions  to  modify  the  puff  and  plume  models  for 
toxic  chemicals  and  the  flash  fire  model  for  flammable  gases  have  been 
made  to  the  VM.  A  general  overview  of  each  modification,  changes  to 
the  user  input,  and  flowcharts  of  major  programs  and  subroutines  are 
given  in  this  chapter. 


PHASE  I 

Executive  Program,  VMEXEC  < Figure  5~1) 

The  main  program  of  the  Vulnerability  Model  (VM) ,  called  VMEXEC, 
initializes  the  data  files  and  controls  the  execution  of  the  various 
simulation  submodels  over  each  time  step  for  each  geographic  cell. 

VMEXEC  is  now  modified  to  call  a  subroutine  DOSAGE  at  the  end  of  the 
program.  The  subroutine  DOSAGE  integrates  the  dose  for  the  puff  and 
plume  models  for  each  geographic  cell.  The  time  step  call  in  VMEXEC 
operates  as  usual  because  models  such  a  pool  burning  still  use  the  time 
step  calculation. 

New  Subroutine,  DOSAGE  (Figure  5-2) 

Subroutine  DOSAGE  computes  the  outdoor  dose  and  indoor  dose  for  each 
cell  for  plume  and  puff  models.  This  subroutine  recalls  some  variables 
from  the  State  file.  The  dispersion  coefficients  for  the  plume  model 
are  from  subroutine  CSSIGS  and  those  for  the  puff  model  are  from  sub¬ 
routine  ISSIGS.  If  the  outdoor  dose  is  less  than  the  threshold  value, 
the  outdoor  and  indoor  doses  are  set  equal  to  zero.  The  indoor  dose  is 
set  equal  to  zero  when  it  is  less  than  the  threshold  value.  This  will 
save  some  computing  time.  Both  outdoor  and  indoor  doses  are  written  on 
a  tape  which  is  transferred  to  Phase  II  for  processing. 

Subroutine ,  FLFIRE  (Figure  5-3) 

The  subroutine  FLFIRE  computes  the  thermal  dose  for  the  burning  and 
cooling  periods  for  personnel  fire  casualties  for  each  cell.  It  also 
computes  the  accumulated  thermal  radiation  for  structural  ignition  at 
each  cell.  Both  thermal  dose  and  accumulated  radiation  for  ignition 
are  written  on  tape  for  Phase  II  processing.  In  addition,  for  the 
secondary  fire  model  the  subroutine  computes  the  distance  from  the  flash 
fire  to  each  secondary  fire  source  and  the  radiation  flux  from  the  flash 
fire  received  at  the  secondary  source  location. 
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Figure  5-1 .  Flow  Chart  of  Executive  Program ,  VMEXEC 
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Figure  5-1  ( continued ) 
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PHASE  II 


For  Phase  II,  no  programming  changes  have  been  made.  Changes  to 
numerical  values  only  have  been  made  to  incorporate  alterations  in 
thermal  criteria  and  probit  equations.  The  flow  charts  for  Phase  II 
are  the  same  as  those  reported  in  reference  [71 . 


USER  INPUT 

Several  additional  user  input  variables  are  necessary  in  modifying 
the  VM.  These  are  given  in  Table  5-1. 


Table  5-1.  ADDITIONAL  USER  INPUT  VARIABLES 


Field 

Number 

Default 

Value 

Unit 

Variable 

Name 

Comment 

2011 

36000. 

erg/cm3 -C 

HEAT 

Summation  of  density 
and  heat  capacity 

2026 

600. 

g/cm3 

RHOT 

Average  density  of 
product 

2033 

100.0 

cm/s 

FLMSPD 

Flame  speed 

1019 

800.0 

°C 

TFLAM 

Flame  temperature 

UPDATE  PROGRAMS 

Two  update  programs  are  given  in  Appendix  B.  The  program  MODVMA 
modifies  Phase  I  of  the  VM,  and  the  program  MODVMB  modifies  Phase  II  of 
the  VM.  The  magnetic  tape,  number  S13346,  contains  this  updated  version 
of  the  VM. 
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Figure  A-2.  GEOGRAPHIC /DEMOGRAPHIC  FILE 
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•Sm  Figures  A-2a  and  A-2b  for  record  formats  for  inter  area  cells  and 
land  area  cells,  respectively. 
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Figure  A-2a 

RECORD  FORMAT  FOR  WATER  AREA  CELLS 


Columns 

Fermat 

Field 

1-8 

(A8) 

‘Cell  Identification 

9-16 

(18) 

‘Latitude,  north 

17-24 

(13) 

‘Longitude,  west  1 

25-28 

(14) 

‘Depth  In  feet 

•Length  in  feet 

35-57 

(23X) 

--not  used — 

58-60 

(A3) 

Direction  of  current  ■ 

61-78 

(18X) 

--not  used — 

79-80  j  (12) 

‘ignition  source 

•Data  used  by  Vulneraoility  Model  for  this  study 


Figure  A-2b 

RECORD  FORMAT  FOR  LAND  AREA  CELLS 


Columns 

Format 

Field 

1-8 

(48) 

•Cell  identification 

9-16 

(18) 

'  -tttude,  north 

17-24 

(18) 

‘Longitude,  west 

25-28 

(14) 

‘Depth  (blank  or  zero) 
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(16) 
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(13) 

Percent  under  18  yr 
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(13) 

Percent  over  62  yr 
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(13) 

‘Percent  sheltered 

44-48 

(15) 

‘Total  housing  units 
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(IS) 

♦Average  value  (in  $100) 

54 

(Al) 

Housing  material 
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(12, Al) 

Number  of  schools 

58-60 

(A3) 

Land  use 

61 

(ID 

Land  uniformity 

62-78 

(17X) 

— not  used — 

79-80 

(12) 

*Data  used  by  Vulnerability  Model  for  this  study. 


_  _  I  — — 


Figure  A-3.  INPUTS  FOR  AMMONIA  SPILL 
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Figure  A-4.  INPUTS  FOR  CHLORINE  SPILL 
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Figure  A-5.  LIST  OF  STATE  FILE  VARIABLES  WITH  CURRENT  DEFAULT  VALUES 


Field 

Number 

Default 

Value 

Unit 

Code 

Display  Name 

Coessent 

1001* •* 

OTW 

01 

CHEN  NAME 

Name  of  cargo 

1002 

200.0 

23 

M0LEC  WEIGHT 

Molecular  weight  of  cargo 

1003 

282. 

06 

BOIL  TEN  LIQ 

Boiling  temperature  of  liquid  cargo 

1004 

0.87 

04 

DENS  LIQ  AHB 

Density  of  liquid  at  ambient  temperature 

100S 

0.261846 

22 

VISCOSITY-BP 

Viscosity  of  liquid  at  boiling  point 

1006 

1.9739 

22 

VISCOSITY-AM 

Viscosity  of  liquid  at  ambient  temperature 

1007 

0.60 

09 

HEAT  CAPC  LQ 

Heat  capacity  of  cargo  in  liquid  phase 

1008 

20.0 

18 

SURF  TENSION 

Surface  tension  of  cargo 

1009 

0.001 

04 

DEN  FUEL  VPR 

Density  of  fuel  vapor  at  boiling  point 

1010 

7.8 

01 

VPE  COEFF  A 

1011 

1443. 

01 

VPE  COEFF  B 

Coefficients  (A,B,C)  of  the  vapor  pressure  equation 

1012 

273. 

01 

VPE  COEFF  C 

1013 

0.1 

09 

HEAT  CAPC  VP 

Heat  capacity  of  cargo  in  vapor  phase 

1014 

136. 

10 

HEAT  OF  VPR 

Heat  of  vaporization 

1015 

4. 

15 

BURNING  RATE 

(Equivalent  to  t5018) 

1016 

800. 

06 

AD  FINE  TEMP 

(Equivalent  to  #5017) 

1017 

0.5 

01 

MOLEC  RATIO 

Not  used  in  VM 

1018 

2.5 

01 

STO  AIR/FUEL 

(Equivalent  to  #5015) 

1019 

800. 

06 

FLAME  TEMP 

Flame  temperature 

1020 

.01 

01 

MOLE  FRACTN 

Water  concentration  below  which  evaporation  is  negligible 

1021 

.87 

04 

LIQ  DENS  BP 

Density  of  liquid  at  boiling  point 

2001 

10000. 

03 

TANK  VOLUME 

Volume  of  tank 

2002 

1000. 

02 

TANK  HEIGHT 

Height  of  tank 

2003 

0.0 

02 

HOLE  HEIGHT 

Height  of  bottom  of  hole  above  bottom  of  tank 

2004 

20.0 

06 

TEMP  START 

Temperature  in  tank  before  discharge 

2005 

0.0 

05 

TANK  PRESS 

Pressure  in  tank  before  discharge 

2006* 

1 

01 

ISO-O,  ADB-1 

Indicator  specifying  isothermal  or  adiabatic  tank  conditions 

2007 

5000. 

08 

INITIAL  MASS 

Initial  mass  of  cargo 

2008 

50. 

02 

HOLE  DIAM 

Average  diameter  of  hole  in  tank 

2009* 

200 

01 

NUM  MASS  INC 

Number  of  increments  used  for  venting  integration  (MOD A) 

2010 

100000. 

02 

OBSVR  DIST 

Distance  from  burning  pool  at  which  flux  is  computed 

2011 

36000. 

27 

HEAT 

Summation  of  density  x  heat  capacity 

2012 

100000. 

02 

COORD  X  «* 

Distance  in  downwind  direction 

2013 

0. 

02 

COORD  Y  •* 

Distance  in  crosswind  direction 

2014 

100. 

02 

COORD  2  •* 

Height  above  ground  surface 

2015 

0. 

02 

HT  OVER  SURF 

Height  of  centerline  of  hole  in  tank  above  water  surface 

2016 

200. 

15 

WIND  SPEED 

Wind  spaed 

2017 

6 

01 

ATMOS  COND 

Atmospheric  stability  flag 

2018 

2 

01 

CHL*1,  RAD-2 

Flag  for  channel  spill  or  radial  spill 

2019 

5000. 

02 

CUR  DIM  POOL 

Current  dimension  of  liquid  pool 

2020' 

50000. 

02 

CHAN'L  WIDTH 

Channel  width 

2021 

0.0001 

16 

HAZARD  CONC 

Not  used  in  VM 

2022* 

0 

01 

FLUX,  CN1,  LM2 

Flag  specifying  conditions  of  heat  transfer  (for  MODD) 

2023 

15.0 

06 

WATER  TEMP 

Water  temperature 

2024 

2.0 

14 

HEAT  FLUX 

Heat  flux  between  water  and  chemical 

*  Denotes  integer  variable. 

•*  Coordinates  of  point  at  which  the  concentration  in  air  is  being  calculated. 


Figure  A-5  (continued) 


Field 

Number 

Default 

Value 

Unit 

Code 

Display  Name 

Comment 

2025* 

1 

01 

CRIT  FLAG 

Flag  for  calculation  of  critical  values  (MOOD) 

2026 

600. 

04 

RHOT 

Average  density  of  product 

2027 

600. 

11 

EVAP  TIME 

Not  used  in  VM 

2028* 

1 

01 

SPILL  INDC 

Flag  indicating  spill  environment 

2029* 

0 

01 

SPL  DUR  INDC 

Indicator  for  duration  of  discharge 

2030 

0.0 

02 

EMMIS  POWER 

Emissivity  power 

2031 

1000000. 

02 

MAX  DIST  HFM 

Not  used  in  VM 

2032 

36000 

11 

MX  TIME  CONC 

Not  used  in  VM 

2033 

100.0 

15 

FLM5PD 

Flame  speed 

2034 

0.0 

02 

MIN  HA2  ZONE 

Not  used  in  VM 

2035 

1000000. 

02 

MAX  HA2  ZONE 

Not  used  in  VM 

2036 

15.0 

06 

TEMP  LIQUID 

Temperature  of  liquid  discharged 

2037 

36000. 

11 

MAX  TIME  CONC 

Not  used  in  VM 

2038 

100000. 

12 

AVG  ESC  RATE 

Average  escape  or  discharge  rate 

2039 

100000. 

02 

CONC  PT  X  ** 

Downstream  coordinate 

2040 

100000. 

02 

CONC  PT  Y  *» 

Cross-stream  coordinate 

2041 

0. 

02 

CONC  PT  Z  »* 

Depth  coordinate 

2042 

600.0 

11 

TIME  CONC  PT 

Not  used  in  VM 

2043 

50. 

19 

DIF  COEF  H20 

Diffusion  coefficient  of  liquid 

2044 

10000. 

02 

RIVER  DEPTH 

Mean  depth  of  flowing  water 

2045 

50000. 

02 

RIVER  WIDTH 

Mean  width  of  flowing  water 

2046 

0.0 

02 

OFF  DIST 

•  Cross-stream  position  of  discharge 

2047 

100.0 

15 

STREAM  VEL 

Velocity  of  flowing  water 

2048 

50.0 

15 

TIDAL  VEL 

Maximum  amplitude  of  tidal  current 

2049 

518400. 

11 

TIDAL  PERIOD 

Tidal  period 

2050 

0. 

u 

PHASE  LAG 

Time  to  next  highwater  slack  tide 

2051 

0.0 

01 

DECAY  COEFF 

Decay  coefficient 

2052 

0.03 

03 

MANNING  FACT 

Manning  roughness  factor 

2053 

0.1 

19 

DIF  COEF  V-A 

Diffusion  coefficient  of  vapor  in  air 

2054 

15.0 

06 

AIR  TEMP 

Air  temperature 

2055 

36000. 

11 

TIME  LIQ  SPR 

Not  used  in  VM 

2056 

36000. 

11 

LIQ  SPR  TIME 

Not  used  in  VM 

2057 

36000. 

11 

TIME  SPL  COND 

Not  used  in  VM 

2058 

90. 

21 

WIND  TOWARD 

Degrees  from  north  toward  which  wind  blows 

3001* 

3 

01 

PROB  TYPE 

Not  used  in  VM 

3002* 

0 

01 

PRINT  FILE 

Logical  unit  for  output  file 

3003* 

0 

01 

PLOT  RAD  FIX 

Not  used  in  VM 

3004* 

0 

01 

NSF 

Secondary  fire  source  indicator 

3005* 

0 

01 

ISF 

Number  of  secondary  fire  sources 

3006* 

0 

01 

ISHLD 

Shielding  situation 

3007* 

0 

01 

NSFIGN 

Number  of  secondary  fire  ignitions 

3008* 

0 

01 

PLOT  CONC  LQ 

Not  used  in  VM 

3009* 

0 

01 

PLOT  SPREAD 

Not  used  in  VM 

3010* 

0 

01 

PLOT  HVP  LQ 

Not  used  in  VM 

*  Denotes  integer  variable. 

**  Coordinates  of  point  at  which  the  concentration  m  water  is  being  calculated. 

A-7 


Figure  A-5  (continued) 


Field 

Number 

Default 

Value 

Unit 

Code 

Display  Nam* 

Comment 

3011* 

0 

01 

OFFER  FILE 

Logical  unit  for  output  file 

4001 

0.0 

08 

tor  MASS  GAS 

Total  anas  of  gaa  which  eacapes 

4002 

0.0 

08 

TOT  MASS  UQ 

Total  msssof  liquid  which  escape* 

4003 

0.0 

03 

TOT  VOL  LIQ 

Total  volume  of  liquid  which  eacapes 

4004 

0.0 

11 

TIME  OF  RFL 

Elapsed  tine  of  release  froa  tank 

400S 

0.0 

06 

MX  TEMP  TANK 

Maximum  temperature  in  tank  during  release 

400b 

0.0 

02 

FLAME  LENGTH 

Length  of  flame 

400? 

0.0 

02 

DLAM  FLAME 

Diameter  of  flame 

4008 

0.0 

13 

FLAME  ANGLE 

Angle  of  flame  from  vertical 

4009 

0.0 

14 

RAD  FLUX 

Thermal  radiation  flux  (pool  burn) 

4010 

0.0 

04 

VAP  CONG 

Vapor  concentration 

4011 

0,0 

02 

1/2  HA2  ZONE 

Not  used  in  VM 

4013 

0.0 

11 

DUR  HA2  CLD 

Not  used  in  VM 

4013 

0.0 

11 

ARRL  TME  HA2 

Not  used  in  VM 

4014* 

0 

01 

IN  OR  OUT 

Not  used  in  VM 

4015 

0.0 

02 

SPILL  S12E 

Not  used  in  VM 

4016 

0.0 

11 

TIME  LQ  EVAP 

Elapsed  time  for  chemical  to  evaporate 

401? 

0.0 

12 

DISSOL  RATE 

Dissolution  rate  of  chemical  on  riverbed 

4018 

0.0 

02 

POOL  FLM  HGT 

Height  of  flame 

4019 

0.0 

03 

VOL  UQ  POOL 

Voluae  of  liquid  remaining  in  pool 

4020 

0.0 

12 

TOT  EVP  RATE 

Average  evaporation  rat*  of  total  liquid  pool 

4021 

0.0 

11 

DISSOL  TIME 

Time  for  all  chemical  to  dissolve  on  riverbed 

4022 

0.0 

04 

LIQ-H20  CONC 

Concentration  of  chemical  in  water 

4023 

0.0 

08 

MASS  VAP  LIB 

Mass  of  vspor  liberated  from  pool 

4024 

0.0 

02 

SAFE  DIST 

Distance  at  which  concentration  1*  less  than  limiting  value 

4025 

0.0 

02 

MAX  DIM  POOL 

Maximum  dimension  of  liquid  pool 

4026 

0.0 

0? 

POOL  SI2E 

Area  of  liquid  pool  on  riverbed 

402? 

0.0 

02 

POOL  LENGTH 

Length  of  pool  on  riverbed 

4028 

0.0 

06 

TCRIT 

Critical  temperature  of  cargo 

4029 

0.0 

18 

SURT 

Intarfaclal  surface  tension 

4030 

0.0 

26 

CSAT 

Solubility  of  substance 

4031 

0.0 

19 

DIFW 

Dlffuaion  coefficient  of  chemical  in  water,  cm*/**c 

4032 

0.0 

11 

SINK  TIME 

Tim*  of  chemical  to  sink  to  bed 

4033 

0.0 

02 

DIST  TRAV 

Distance  traveled  by  chemical  to  reach  bed 

5001* 

0 

01 

MESSAGE  FLAC 

Flag  vhlch  control*  output  Mssages 

5002* 

0 

01 

MISCIBLE  XND 

Miscibility  indicator 

5003* 

0 

01 

REACTIVE  IND 

Reactivity  indicator 

5004* 

0 

01 

TOX.ASPH  IND 

Toxicity  and .asphyxiation  indicator 

5005* 

0 

01 

LIQ  CONC  IND 

Liquid  toxicity  indicator 

5006* 

0 

01 

IGNITION  IND 

Flag  which  indicates  type  of  ignition 

500?* 

0 

01 

IGNITN  CELL 

Cell  in  which  ignition  first  occurred 

5008* 

0 

01 

ICNITN  CODE 

Ignition  code  of  the  ignition  cell 

5009* 

0 

01 

SPILL  CELL 

Cell  in  which  discharge  occurred 

5010* 

0 

01 

FLAG, ISO, CS1 

Flag  for  use  of  puff  or  plum*  equation 

‘Denote*  Integer  variable. 
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Figure  A-5  ( continued ) 


Field 

Number 

Default 

Value 

Unit 

Co^ie 

Display  Name 

Comment 

5011 

1.0 

0*1 

DENS  WATER 

Density  of  water 

5012 

-40. 

00 

FREEZING  PT 

Freezing  point  of  chemical 

5013 

0.0 

01 

CONCi  DO  LIM 

Lower  flammability  limit,  percent 

5014 

0.0 

01 

case,  UP  LIM 

Upper  flammability  limit,  percent 

5015 

0.0 

01 

AIR/FUEL  RAT 

Stoichiometric  air/fuel  ratio 

5016 

0.0 

10 

HEAT  COMBUST 

Heat  of  combustion 

5017 

0.0 

06 

ADFLAME  TEMP 

Adiabatic  flame  temperature 

5018 

0.0 

IS 

BURNING  RATE 

Burning  rate  of  chemical 

5010 

0.0 

01 

MOLES  OXYGEN 

Moles  of  oxygen  per  mole  of  fuel 

5020 

282. 

06 

FLASHPOINT 

Flashpoint  of  chemical 

5021 

1.4 

01 

SPC  HEAT  RAT 

Ratio  of  specific  heats 

5022 

1013250. 

05 

PRESSURE  AMB 

Ambient  atmospheric  pressure 

5023 

C.O 

04 

CONC.  LO  LIM 

Lower  limit,  flammable  vapor  concentration 

5024 

0.0 

04 

CONC,  UP  LIM 

Upper  limit,  flammable  vapor  concentration 

5025 

0.0 

04 

CONCi  STOICH 

Stoichiometric  vapor  concentration 

5026 

0.0 

05 

VP  PR  TANK 

Vapor  pressure  in  the  tank 

5027 

0.0 

08 

LIQ  REM  TANK 

Total  mass  of  liquid  remaining  in  tank 

5028 

0.0 

02 

SIGY , DISPERS 

Horizontal  dispersion  coefficient  (MODC) 

5029 

0.0 

02 

SIGZiDISPERS 

Vertical  dispersion  coefficient 

5030 

2.75 

01 

EXP.  TOX  SUM 

Exponent  in  weighted  sum  for  toxicity 

5031 

5032 

0.0 

0.0 

01 

01 

Tl.  COEFF  A 
Tl,  COEFF  B 

Coefficients  of  probit  equation  Tl 

5033 

5034 

0.0 

0.0 

01 

01 

T2,  COEFF  A 
T2,  COEFF  B 

Coefficients  of  probit  equation  T2 

5035 

0.0 

01 

T3.CONC  PPM 

Irritation  threshold  for  T3 

5036 

100.0 

01 

T4.  COEFF  1NG 

Coefficient  of  ingestion  T4 

5037 

0.0 

01 

CUAL  STD  AIR 

Not  used  in  VM 

5038 

0.0 

01 

FRAC  POP  SHL 

Fraction  of  population  sheltered 

5039 

0.0 

02 

FFDIA 

Flash  fire  diameter 

5040 

0.0 

11 

TIME  AT  IGN 

Time  at  which  ignition  occurred 

5041 

0.0 

02 

PUFF  CENTER 

Downwind  location  at  center  of  vapor  puff 

5042 

0.0 

14 

FL  FIRE  RADN 

Intensity  of  flash  fire  radiation 

5043 

0.0 

11 

FL  FIRE  TIME 

Effective  duration  of  flash  fire 

5044 

0.0 

08 

MASS  VAP  EXP 

Mass  of  cargo  vaper  which  exploded 

5045 

0.0 

01 

DCT  VAP  EXP 

Mass  exploded  given  as  percent  of  total  mass 

5046 

0.0 

20 

EXPL  YIELD 

Yield  of  the  explosion 

5047 

0.0 

2S 

TNT  EQUIV 

Yield  equivalent  in  short  tons  of  TNT 

5048 

0.0 

U 

PL  BURN  TIME 

Effective  duration  of  pool  burn 

5049 

0.0 

01 

RAD'  COORD 

Radius  of  vapor  cloud 

5050 

0.0 

11 

CURRENT  TIME 

Current  value,  elapsed  time 

6001* 

0 

01 

TIME, SEC, BEG 

Time  in  seconds,  begin  loop  1 

6002* 

0 

01 

TIME, SEC, END 

•  Time  in  seconds,  end  loop  1 

6003* 

0 

01 

TIME, SEC, INC 

Time  in  seconds,  increment  value  for  loop  1 

‘Denotes  integer  variable. 
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Figure  A-5  ( continued ) 


Field 

Number 

Default 

value 

Unit 

Code 

Display  Name 

Comment 

6004* 

0 

01 

TIKE, MIN, BEG 

Ties  In  minutes,  begin  loop  2 

6005* 

0 

01 

TIME, MIN, END 

Time  in  minutes,  end  loop  2 

6006* 

0 

01 

TIKE, KIN,  INC 

Time  In  minutes,  increment  value  for  loop  2 

6007* 

0 

01 

TIME, KIN, BEG 

Time  in  minutes,  begin  loop  3 

6008* 

0 

01 

TIME, MIN, END 

Time  in  minutes,  end  loop  3 

6009* 

0 

01 

TIME, MIN, INC 

Time  in  minutes,  Increment  value  for  loop  3 

6010* 

0 

01 

SPILE  LAT 

Latitude,  north,  of  spill  site 

6011* 

0 

01 

SPILL  LONG 

Longitude,  west,  of  spill  site 

‘Denotes  integer  variable. 

Figure  A-6 

INPUT  DATA  FOR  LNG  SPILL 
(FLASH  FIRE  TEST  RUNS) 


•1001  two - 

1002  16.04 

2011 — -.13726*05 - - - 

2026  .12346-02 

2001  — -.30006*11 .  - . 

2002  1500. 

-2003 - 100*0 - 

2004  -151.0 

200S--  .20276 *07 - - 

2006  1.000 

2007  - .12506*11 -  - - 

2008  100.0 

-2014 - HJOtO - 

2015  0. 

2016  —  400. 0-- . . . 

2017  4.000 

2019”  -2. 000 - - —  -  • 

2020  .6096^*05 

-202-1 - r 1-2466-0  3 - 

2022  1.000 

2023  —  20.00 .  . - 

2025  1.000 

2027—  180.0 - - - 

2033  500.0 

-2054"- — 15t00  — . . . . . . . . . 

2058  30.00 

-3004  1.000—  -  • 

3006  2.000 

5001  •  0. -  •  •• 

5002  0. 

■5003 - WOOO - 

5004  0. 

5010  — 1.000  -  • 

5011  1.000 

5012  -182.0  -  -  - 

5013  5.400 

•5014 - !<r*00 - 

5016  -.13006*05 

5017  •  1100. - - 

5018  .20836-01 

5019  2.000 

5020  -161.0 

■5021 - 1-.400 - 

5022  ,10136*07 

5036  -  0.  ----- 

5037  0. 

5038  .5000 
6001  0. 

6002 - 0* - 

6003  0. 

6004  2.000  - 

6005  10.00 

6006  2.000  -  - 

6007  0. 

6008  0. 

6009  0. 

6010  .29546*06 

6011  .90076*06 

0  0.  -  - 


UPDATE  PROGRAM#  MODVMA 


•IDENT  VMEXtO 
•DELETE. VMEX1. 2 

1TAPE9. TAP2l0.TAPEl2*TAPEl3,TAPEl4.TAPEiS*TAPt24.TAPt33. 
1TAPE34»TAPE30«SECFHE) 

•1NSLRT  .VMEXEC.  16 

COMMON/PAOOS/Cl  ( 400 )  .02(400) 

•lNSEHT«VMExEC.144 

IFUTOX  .EQ.  1)00  TC  901 
•INSERT. VMEXl .48 

WHITE  <b»704) 

704  FORMAT  (lhl»10X»7HCELL  NO. 8X . 1 1 hPEOPLE  L'OSt  ♦  1  OX  ♦  1  OHhOuSE  DOSE) 

00  70  l*l.NCELL 

IF(GUI)  .EQ.  0.0  .AND.  02(1)  .EQ.  0.0)00  TO  70 
WHITER. 702)I.Ul(I).Q2(I> 

70  CONTINUE 

702  FORMAT (10X. I5.5X. 1E16.4.5*. It  lb. 4) 

•0ELETE.VMEXS.36 
C  WR I TE ( 6 .91 94 ) 

•UELETE.VMEX1.I44 

C  WRITE  (6.9195JLC.  ( ( IP.TINCUP)  «TFRI  (1P.LC)  >  v*iP*l»NlNC 
•DELETE .VMEXl.l4ri.VMtX 1.152 

C9194  h OHmAT ( » i  riMt  INCREMENTED  RADIATION  FLUX  AT  EACH  CELL*.//. 

C  1  «  CELL  iNCKEMtNT  tFFECTIVE  CURA1 ION  HAU1 AT ION* ./ . 

C  2  *  NO.  NO.  (SEC)  FLLA ( J/M2/S) * ) 

C919S  FORMAT ( (2X.13.SX.13.11A.F8.2.7X.E12.H.21 (/♦ 10X. 13. 1 lX.Fo.2. 

C  1  7X.E12.4)) 

•INSEKT.VMEXEC.2J1 

CALL  FKCL(4016*TEVAP. Ib.lH) 

CALL  FPCL (4023.TMVAP. Ib. 1H) 

CALL  FHCL (20lto.UWlND.is.IW) 

CALL  FHCL (4030. TREND. IS. 1R) 

CALL  1HCL (5004.ITOX.IS.IK) 

IFUTOX  .EQ.  1  ) CALL  UUSAOt 

•IDENT  ERFFA 
•OELETE.EHF .22 

IF ( V  .GT.  6.0  .OR.  Y  .EQ.  6.0)tX=0.0 
IF ( Y  .LT.  6.0)EX*tXP(-Y*Y) 

•IDENT  DOSAA 
•INSERT  OMSOEG.ll 

SUBROUTINE  DOSmUE 

COMMON/CELL/IXLAT(400)  •,  1YLON(400)  .CuEP(AOO)  .CLEN(AOO)  .10COOE(400) 
lX(40Q)*Y(400).  IGCtLL  (*.00 )  » NCtLL  ♦  1  SPLAT  .  1 SPL  ON.  WINDED 
DIMENSION  CX(SOO) .DX(SOO) .EX (SOU) 

DIMENSION  005(400). OOSlN(400) 

C  DOS  *  OUTDOOR  DOSAGE 

C  D0SIN  *  INDOOR  DOSAGE 

CALL  FHCL (1 002 .WTMOL.  IS.  IH) 


CALL  FRCU2014.Z.1S.IR) 

CALL  FRCL <20l6.UWlNC.Ib.lR) 

CALL  FRCLt2019.SiZE.lS.IR) 

CALL  FRCU2020.CHNLWtib.lR) 

CALL  FRCL (4001.TMGtIS.lR) 

CALL  FRCLt4QIfc.TtVAP.Ia.IR) 

CALL  IRCU2017.1AC.IS.1R) 

CALL  IRCLt201fl.10lM.IS.lfi) 

CALL  FRCLt205A.AlRTfM.iS.iR) 

CALL  iRCLtSOlO.  SCS.IS.1R) 

CALL  FRCLt5030.CP0w.lS.iR) 

CALL  FRCL (5031.TlA.IS.iR) 

CALL  FRCU4023.TMVAP.Ia.lR) 

CALL  FRCLt5032.TlB.lS»lR) 

CALL  FRCL (S033.T2A. IS. 1R) 

CALL  FRCU5034.T2B.IS.1R) 

CALL  FRCLt5035.T3A.IS.iR) 

OATA  PI/3.141SV/ 

ThRS  =  3.0 
NCOD  =  1 
RAU  s  SIZE 

THRS  s  THRESHOLD  FOR  IRRITATION.  PPM 
SQ2  =  SORT (2.0) 

ALF  =  28.966»l.E6/(.OOi225«w':MOL) 

ALF  IS  A  CONVtxSlON  FACTOR  CHANGlNti  G/CC  TO  PPM  IN  VOLUME 
28.966  IS  AIR  MOL  wT.  .001223  IS  AIk  DtNSHY  AT  SAMt  CONDITION 
AS  THE  ChtMICAL 
VWlND*UwlNO<*.  022369 

VWlNO  IN  MPH 
TEMU=AdS(25.0-«IRTEH) 

R=< .25*.02165«VWlNO*.0uB33*TtM0)/60.0 
R  IS  THt  INFILTRATION  PER  MINUTE 
IF  (TMG  .EQ.  0.0)TMG*TMVAP 
IF ( 1 0 1 M  ,EQ.  1  ) RAQ=SOkT  t S I Zt#CHNLW/PI ) 

IFtlSCS  *EO.  0 ) GU  TO  lul 
PLUME  MODEL 

OB  a  TMG/(PI*TtVAP*UwIiNO) 

DOS  IS  IN  PPM  ANO  MINUTE 
DO  5  Is! .NCELL 
OOS(I)=0. 

OOS IN  1 1 ) =C0S  ( I ) 

IF  < X ( I )  .LT.  0 • ) GO  TO  b 

XA  *  XU)  ♦  10.0«RA0 

CALL  CSSIGStXA.IAC.SIGY.SIGZ) 

IF  tSlGZ  .GT.  100000. ) S IGZ  =  10O000 . 

CC  IS  THE  MAXIMUM  VAPuR  CONCENTRATION 
CC  s  ALF*UB/tSIGY*SlGZ) 

CEs0.5#((Y(I)/SIGY)**2*<Z/SloZ)##2) 

CAX»CC**CPOW 

CAYsCPOw«CE 


B-2 


CG=0. 

IF (CAY  .LT.  23.0)  CG=CAX/EXP(CAY) 
OOS(I)=sCG*TeVAH/60.0 

IF ( T3A  «NE «  0.0  .ANG.  CC  .LE.  T3A) DOS ( I ) =0.0 
IF ( T3A  .EG.  0.0  .anC.  CC  .LE.  ThRS) DOS ( I ) =0 .0 
CALCULATE  PLUME  INDOOR  DOSAGE 
OUSIN(I)  =  DOS ( I ) 

IF (DOS ( I )  .EQ.  0.0)GO  TO  5 

cxu)  =cpow 
OX(l)  =  l.o 
CRT  =  .0001*CPOW 
YE  =  EXP(-R*TEVAP) 

AXA  =  CPOwMl.O  -  YE) 

00  91  J=1.100 
XI  =  J 

DX(JM)  =  <XW.0)<MXIM.0)*DX<J) 

CX  <J*1  )=CX  ( J)MCPOw-XI) 

AXB=(-l)**J*<CX(J*l)/DX(o*l) )*(l.-Yt**(J*l) ) 
AXA=AXA+AXB 

IF1AXB  .LT.  CRT)G0  TO  93 
91  CONTINUE 
93  CONTINUE 

DOS  IN ( I ) =CG* ( TEVAP-AXA/k) /bO . 0 
5  CONTINUE 
GO  TO  900 
101  CONTINUE 

CALCULATE  PUFF  COSAGE 

BA=^.0*TMG*ALF/(2.0«PI)*»1  .5 

BC  IS  ThE  MAXIMUM  CONCENTRATION  IN  THE  CELL 

BB  =  (2.0/CPOW)**0.5 

DO  15  1=1 .NCELL 

DOS ( I ) =0 . 

uOSIN ( I ) =0 . 

IF ( X  ( I )  .LT.  0.)GO  TO  IS 

XA=X(I)M0.0«HAD 

CALL  iSSIGStXA* IAC»SIGY.SIGZ) 

SIGX  =  SIGY 

8C=(BA/(SIGX*SIGY*SIGZ) ) 

BD=BC**CPOW*BB*S IGX/Uto INC 

BE=CP0W«(  (Y(I)/SIGY )*«<:♦  (Z/S10Z)#*2>/2.0 

IF (BE  .GT.  23. ) GO  TO  lb 

IF(BE  .LT.  23.0  .ANC.  oE  .GT.  1  .E-2B) 6F  =BD/EXP (BE) 
IF  (BE  .LT.  1.0E-2b)bF=oD 
BG= (CPQ*/2 .0)**0.S*X(I) 

Bh=(PI**0.5/2.0)*(1.0*EkF (BG/SIGX) ) 

DOS ( I ) =BF *8P/G0 • 0 

IF ( T3A  .NE.  0.0  .ANC.  dC  .LE.  T3A) DOS ( I ) =0 .0 
IF  (  T3A  .EG.  0.0  .AND.  dC  .LE.  3. 0) DOS ( I ) =0 .0 
TO  CALCULATE  INDOOR  COsAGE 
IF  (DOS ( I )  .EQ.  0.01G0  10  lb 
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EA=(SQWT(2.*PI )*SIGX»ft/(2.*60.*UWlN0> )**CPOw 

£8  =  CP0V**  (X  ( I  >*k/(fcO.*U*INO>*(SlGX*P/(6Q.*U#lND>  )**2/2. 

EX=£B-8E 

IF (ASS (EX)  .Lt.  l.E-25)ECsBCtt*CP0W*tA 

IF  <  ABS ( EX )  .GT.  i.E-25) tC=bC*«CPOW*£A*EXP<tA) 

E0=(X(i)*(W*SioX**2)/<b0.*U*lN0>  >/(SG2<»SIGA) 

IF  <EQ  .LT.  5.0) £F=ERF (tD) 

IF  <£0  .GE.  5 • 0 ) EF  =  1 .0 

FA=(X(I)  ♦(SIGX**2)*F/(bO.«U*iNC)  )/U*lNC 
I r2=F  A 

FB=F  A-5.0«S()2*bI6x/trflND 
ITl=Fd 

FC=FA*5. 0*502*5 IGX/LwInO 

IT3“FC 

XT=ITI 

SUM-0 • 

DO  20  i\aITl.IT3 
EI=EXP(-CPOW«k*XT/60.U) 

EM=EXP  (~CPO**ttP*  (XT  *1 .  )  /  fc 0 . 0  ) 

IF (XT  .GE.  I T2 ) GO  TC  2b 
EG=ED-U*I.\0*XT/  (S024»SIuX) 

IF (EG  .LT.  5.)th=tPF(E0) 

IF (EG  .GE.  5.)tH=l. 

EJ=EI'MEF-Eh)  o^CPO* 

XT=XT*1. 

EK=ED-U*lND*XT/(SQ2«SIoX> 

IF (EK  .LT.  5.)tL=EPF(t»0 
IF ( £K  .GE.  5.)tL=l. 

EN=EM*(EF-EL)«*CPO* 

SUM=SUM*(EJ*EN)*l./2. 

GO  TO  20 
25  CONTINUE 

EG=UWlNC*XT/<SU2«SIGX)  -  ED 
IF  (EG  .LT.  5.  )  th=EfcF  (EG) 

IF (EG  .GE.  5.)EH=1. 

EJsE  I*  (EF*EH)  <>*CPOW 
XT=XT* 1 . 

EK=UWlND*XT ',Su2*SIGX)  -ED 
IF  (EK  .LT.  5.)tL=fcHF(EK) 

IF  (EK  .C-E.  5.)EL  =  1. 

ENsEM*(EF*EL>**CPOw 

SUM=SUMMEJ*EN)*I./2. 

20  CONTINUE 

EP=CPOW*P«(FA>5.0*S(J2*SlGx/UtNiND)/60.0 

EQs2.0**CPOW»EP*b0.0/<CPCw*H) 

DOSIN(I)=EC*(SUM+EG) 

15  CONTINUE 
900  CONTINUE 
REWIND  33 
WHITE (6*7A3) 
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743  FORMAT ( IH1 . 10X*4HC£LL.tiX. 12HU0SE  0UTSICt.6X. 1 lhDOSt  INSIDE/) 
DO  752  1  =  1.  NCELL 
WRITE(33»741)CC3(I) »OOS IN ( I ) 

752  CONTINUE 

DO  757  1  =  1  .NCELL 

IF ( DOS ( I )  .tO.  0  • )  GO  TO  757 

WRITE (6.742) ItUOS(l) .DUSlN(I) 

757  CONTINUE 

742  FORMAT (1  OX. l4.dX*El2.4»5X«E12.4) 

741  FORMAT (2E16. 6) 

RETURN 

ENO 

*  I  DENT  MOOZA 

•Insert. mccz2. i 

COMMCN/RACOS/Gl (400) .02 (400) 

•ident  flaba 

•INSERT. FL^IkI. 4 

COMMON/RAOOS/O I (400) .02 (400) 

•DELETE. FLFIRl.StFLFlRE.il 
DATA  PI/3. 1415*2567/ 

•DELETE. FLF IR1. 6. FLFIR 1.30 

DATA  SIGMA /5.6o*7t-5/«MLF/0.*/ .RAM/ 1 • 34 t ♦ 7/ 

C  RAM  IS  ThE  CRITERION  FuR  WOOD  BRUN.  tRfc/CM2-b 
CALL  FRCL(20ll.htAT.IS.I«) 

CALL  FRCL(2026tRHOT.IS.IR) 

CALL  FRCL(1019.TFLAMtIStlR) 

CALL  FR CL (5023 .CONCLO. IS. 1R) 

CALL  FRCL(5025tR0F.IS»lR) 

CALL  FkCL  (2033.FLMSPD. IS. IR) 

C  Ht AT  IS  THE  SUM  OF  CENsITY  X  CP 

TF=TFLAM 

C  RHOT  IS  THE  SUM  OF  PROOUCT  DtNSITY 

ABAsTMV/ ( (2.0*Pl ) **1 •5*SIGX*SlGY*SlGZ*KL) 

C  KL  IS  THE  lower  flammable  LIMIT 

RL=(3.O«TMV/(4.0*PI*ROF> )*«(l./3.) 

C  ROF  IS  THE  FUEL  DtNSITY 

QIR  =  SIGMA*(TF«*4  -  TA**<t)  •  TF 
C  SIGMA  IS  STEFAN-UOLTZMANN  CONSTANT 

C  ChANGE  HC0M8  FROM  CAL/G  TO  ERG/G 
HBURN=HCOMd*4.184E»07 
ACB=-ROF»TA*HBURN-(TF-TA)*hEAT«TA 
HCOMB  IS  IN  NEGATIVE  VALUE 
SV=GIR/ACB 

SV  IS  THE  FLAME  VELOCITY 
FLMSPO  IS  THE  INPUT  FLaME  SPEED 
TOTMA  IS  THE  AIR/F'UEL  MASS 
T0TMA=TMV*18.1B69 

R6=(3.0*TMV*TF/(4.0*PI •RCF  *T  A ) )**(l./3.) 

C  MB  IS  THE  MAX  FIREBALL  RADIUS 
SV=FLMSP0 


TFOUR=«e.SV 

CALL  FSV  (5043*  TFOURfA) 

X=XC<IGMU) 

Y=YC (  IONIC) 

C  IGNID  IS  THE  NUMBER  OF  IGMTtO  CELL 

TR  =  T  A/TF 

C  HERE  TA  AND  TF  ARE  IN  rvtLV  IN*  SEE  MOUZ  SUBROUTINE 
CD=.02066b»TA<»*  (23./9.) 

Ct=9.U»  <ALF**9*SI0MA)**  < 1 ./3.>*HEAT«TA 

CF=CE*R8/(TF»<Hll,/9.)  ) 

BC=ALF«SIGMA*TM*a 
bA=bC«*  <4./3. ) 

EA=TA*RB<*REAT/<3.*SIGMm*TF*<M4./3.)  > 

00  5  1  =  l  .NCELL 

01  (I)=0. 

02 ( l ) =0  * 

CONC=QDATAU> 

OIS=(  (XC(l)-X)**2>YC<n**2)<>*0.5 
IF  (I  .NE.  IGNlU) GO  TO  ibO 
WRITE (6t9l ) I *CONC 

91  FORMAT  ( *0  CELL*tI4.<*  IN  FLASH  FIRt  ♦  CONC=  '*♦£12. 3) 

150  CONTINUE 

dB=<Rb/CIS>tt*(b./3.) 

IF (RH  .GE.  0 IS ) Bb  =  l  . 0 
QA=3.0*BA«*Bb*Rb/ ( 1 1 .*FlMSPD) 

C  CHANGE  I  TO  J0ULE/M2 ♦  SEE  ORIGINAL  tUUATION 

CA  =  l  ./23.-TR*«A/  13. -Tft»*d/49.-TR*<*  12/05. 
Qb=CF*6B*(CA*CU> 

QlCDsClA 

C  0  IS  THE  TOTAL  RAOIATIUN  IN  CELL  1 
I F  ( Q 1 ( I  >  .LT.  i  .E*3)G0  TO  5 
C  CALCULATE  STRUCTURE  RAC1ATI0N  HEAT 
IF  1 1  .NE.  IGNIOGO  TO  ob 
02 ( I ) =2  »0£ ♦  13 
GO  TO  5 
66  CONTINUE 

RlM=(RAM*ClS«<*2/oC)«*0.5 
C  HIM  IS  ThE  CRITICAL  KAOICS  FUR  WOOD  BURN 
OC  =  0. 

DA=(BC/0IS»*2>*M  .5/<A.*SV) 

0B=RB«*2-RIM**2 

IF (DB  .LT.  0  • )  GO  TO  56 

0C=RB*D8**1  .5-1 .5«HlM**2«Rb*UB**0.5 

DE=0C*1.5<*RIM#*4*AL0G< (RB ♦Ub**.5 ) /RIM ) 

QC=OA*GE 
56  CONTINUE 
00  =  0. 

E1=RAM*CIS*«2*TF** (2./J.) 
T2=(El/(ALF*SlGMA*Re**2>)**(3./14.) 

E J=T  2/TF 
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IF  (EJ  .GT.  1 *0 ) GO  TO  5 
EK  =  T  A/T  2 

EL=0.1*3.*EJ*«  < 14./3.)/8.  -E J** (26./3. > /48. 
EM=EL-0.5*£K**4O.*EJ««<  14./J.)«EK**4 
EN=.4556-.4i88*EK»«4-.U5358'»tK*«a 
EP=3.*(TF*«(10./3.)*EM-T2*»(10./3.)*EN> 

EQ= ALF ttS  IGMA*kb#*2/ (0 1 b*w2 ) 

QO=tQttE Ao£P 
U2<I>=UC*Q0 

C  02  IS  THE  STRUCTURE  EIuIATlON 
5  CONTINUE 
REWIND  34 
00  10  I=1.NCELL 
WRITE<34»70i>Ul<n.C2U> 

10  CONTINUE 
701  FORMAT (2E16. 5) 

•DELETE .FLFIR4. 2 

IF(Q2(lC)  .GT.  0.0) NPF LAG= 1 
•DELETE. FLFIR5. IS. FLH IR5.25 

IF  (C  1ST  .GE.  RL ) GO  TO  bl 
QSR  <  IF ) =G IR 
GO  TO  55 
51  CONTINUE 

0SR(IF)=<RL««2/(2.*DIST**2> >*UlR 
•IDENT  STRUC 
•INSERT. FLFIR1. 44 

COMMON/RdDOS/Gi (MOO ) .02 (400 ) 

•INSERT. FLFIR1. 57 

CALL  IkCL(5006.IIGN,IS.IR) 

IF ( 1 1GN  .GT.  1 ) GO  TO  94 1 
F  =  0. 

FLIT=0. 

I F (G2 (LC )  .LT.  1.3fc4E*12)RETUNN 

Fsl.O 

FL  IT  =  1 .0 

RETURN 

941  CONTINUE 


UPDATE  PROGRAM,  MODVMB 


*  IDENT  PhASlA 
•DELETE.PH2.1 

lTAPE23«TAPE24»TAPE25tTAPEa6*TAPE33tTAPE34»TAPE6=OUTPUT) 

«INSERT*PHASEII.9 

COMMON/ TOX/UOS  (400)  *D0SIN(400) 

COMMON/HAGOS/U 1 (400) *02(400) 

•insert .phase 1 1 .60 

IF ( ITOX  .NE.  1>G0  TO  7/a 
REWIND  33 
DO  1  1  =  1 *NCELL 

READ (33*777)  DOS ( 1 ) «DOSlN ( I ) 

IF (EOF ( 33 ) ) 965 » l 
1  CONTINUE 
777  FORMAT (2E16.6) 

776  CONTINUE 
REWIND  34 
DO  S  1=1 , NCELL 
RE AD( 34*740)01(1)* 02 (I) 

IF (EOF (34) )965»5 
S  CONTINUE 
965  CONTINUE 
740  FORMAT (2E16. 5) 

*UELETE*PHASEII.74 

130  IF ( ITOX  .NE.  1 ) CALL  PRLONC ( W TMOL *COPW) 

»0ELETE*PhASEII.96 

*  10HF IHtbALL  ./A.12HP00L  BURNING) 

*  INSERT .PHASE  11.141 

IF ( ITOX  .EQ.  1)G0  TO  275 
•IOENT  SADTOX 

•DELETE  SA0TA1.5*SACTA1.35 

COMMON/ TOX/DOS (400) *DOSIN(400) 

PCT1  =  0. 

PCT2=0  • 

PCT3=0. 

IF(DOS(I) .UT.O.O.OR.DOblN(l) .LT .0. ) WRITE (6*431 ) I *CCS( 1 ) .DOSIN ( I ) 
IF (DOS ( 1 ) .LT .0 .0 .OR.DOSIN ( I ) *LT .0. ) 00  TU  900 
431  FORMAT (I1Q.2E16.4) 

IF (DOS ( I )  .EO.  0.)G0  TU  40 
PRT1  =  T1A*T1B«AL0G(D0S(1)  ) 

PCTl=XNORMA(PKI  1-5.) 

PRT2  =  T2A*T2B*ALOG (DOS ( 1 )  ) 

PCT2=XN0RMA (PHT2-5. ) 

PCT2=PCT2-PCT1 

IF (PCT2  .LT.  0 . ) PCT2=0 • 

IF (DOS ( I )  .LE.  T3A ) GO  TO  40 
PCT3=1 •  -PCT1-PCI2 
.  40  PARR ( 1 ,KQUT) *PCT  1 
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PARft(2*K0UT)=PCT2 
PARR (3»KCUT)=PCT3 
C  EVALUATE  INCOCk  deaths 

PCT 1 8=0 • 

PCT19=G. 

PCT20=0 • 

IF(DOSIMI)  .Eu «  0.)G0  TO  75 
PRT18  =  0  • 

PRT 19=0  • 

PRT20=0 • 

PRT18=TlA*Tlti*ALCG(C0SIN( 1 ) ) 

PCT18  =  XN0RMA (PHT  18-5. ) 

PRT19  =  T2A*T2titfAL0G(C05ifK  I)  > 

PCT19  =  XNCRMA(PRT  19-5.) 

PCT19=PCT19-PCTle 

IF (PCT 19  .LT.  0 •  ) PCT 19=0 . 

IF (OOSIN ( I )  .Lt.  T3A ) GO  TO  75 
PCT20=i .-PCT18-PCT19 
75  PAKR(ia.K0UT)=PCT18 
PARR(19,(<CUT)=PCT19 
PARR  (20 .rtOUT )  =PCT20 
900  CONTINUE 
* IOENT  SADFA 
•INSERT. SA0F2. 5 

COMMON/R AOOS/0  1  ( 0 0  )  .02  (400) 
•DELETE »SAQF2.0»SADF£.  13 
PCF 1=0. 

PCF2  =  0  • 

Q=Q1(I)/1.0E*08 
IF (Q  .ED.  0 • ) GO  TO  200 
•DELETE.SADF2.14 

PRF1=-12.8»2  »5to*AL0G (0) 
•DELETE.SA0FA.4 

PRF2=-39.83O.Ul86«AL0G(G«l  .uE^4) 

•INSERT.SA0FA.7 
200  CONTINUE 
•10ENT  SADSS 
•insert. saos2. 6 

COMMON/ R ADOS /Ul  UUO)  .02  (400) 
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METRIC  CONVERSION  FACTORS 


